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Prefacio 


Hace unos años, la Universidad Social me invitó a dictar 
dos conferencias sobre el método de Montecarlo para los 
oyentes de la Facultad de Técnica de cálculo. Yo acepté, 
Poco antes de la primera conferencia me enteré con espanto 
de que los oyentes, en su mayor parte, ignoraban la Teoría 
de las probabilidades. No había modo de retroceder. Sobre 
la marcha hubo que incluir en la conferencia un apartado 
adicional con el fin de exponer a los oyentes las nociones 
principales de la Teoría de las probabilidades. 

El contenido de estas conferencias, dictadas a lo largo 
de varios años, con el tiempo iba estabilizándose. En la 
edición presente se ha conservado también el apartado adi- 
cional ($ 2) a] que quiero dedicar unas palabras. 

Cualquier persona habrá empleado alguna vez las pala- 
bras «probabilidad» y «variable aleatoria». La idea intuitiva 
de la probabilidad (en tanto que frecuencia) corresponde 
más o menos al contenido verdadero de este concepto. En 
cambio, la idea intuitiva de la variable aleatoria está muy 
lejos, como regla, de la definición matemática. Por eso, en 
el $ 2 damos por conocido el concepto de probabilidad y 
explicamos sólo el concepto más complejo de variable alea- 
toria. Este parágrafo no puede reemplazar el curso de la 
Tooría de las probabilidades: la exposición es simplificada 
y ho trae demostraciones. Sin embargo, permite hacerse 
una idea sobre las variables aleatorias suficiente para 
comprender los aspectos elementales del método de Monte- 
carlo. 

El objetivo principal de este libro es sugerir a los espe- 
cialistas de las más diversas ramas la idea de que en el 
campo de sus actividados existen problemas que se pueden 
resolver por el método de Montecarlo, 

Los problemas que se analizan en las conferencias son 
variados y suficientemente simples. Sin embargo, no pueden 
abarcar, claro está, todas las esferas de aplicación del mé- 


todo. Me limitaré a dar un ejemplo. En el libró no se habla 
para nada de la medicina. No obstante, los métodos del $ 7 
permiten evaluar las dosis en la radioterapía. Teniendo un 
programa de cálculo de los rayos absorbidos por distintos 
tejidos del cuerpo, se puede dosificar y orientar la radiación 
del modo más eficiente cuidando de no dafiar los tejidos 
sanos. 

Este libro abarca todo cuanto se irató en las conferen- 
cias; contiene una exposición más detallada de los ejemplos; 
además, se han agregado los parágrafos 9 y 10. 


1, Sóbol 


INTRODUCCIÓN 


$1 
GENERALIDADES 
DEL MÉTODO 


El método de Montecarlo es un método numérico que 
permite resolver problemas matemáticos mediante la simu- 
lación de variables aleatorias. 

1.1. Orígenes del método de Montecarlo. Se considera 
como fecha de nacimiento del método de Montecarlo el año 
1949 en el que apareció el artículo titulado «The Monte Carlo 
method»). La creación de este método suele ligarse a los 
nombres de los matemáticos norteamericanos J. von Neu- 
mann y S. Ulam. En la Unión Soviética los primeros artícu- 
los dedicados al método de Montercarlo aparecieron en 
1955 y 1956 2). 

Es curioso que la base teórica del método era bien conoci- 
da desde hace mucho tiempo. Es más, algunos problemas de 
la Estadística se resolvían a veces empleando las muestras 
aleatorias, o sea, aplicando de hecho el método de Montecarlo. 
Sin embargo, hasta la aparición de las máquinas calculadoras 
electrónicas (MCE), este método no encontraba aplicaciones 
suficientemente amplias ya que la simulación a mano de 
variables aleatorias constituye un proceso muy laborioso. 
Es decir, la aparición del método de Montecarlo en tanto 
que un método numérico de gran universalidad se hizo po- 
sible sólo gracias a la creación de las MCE. 

El nombre de «Montecarlo» se debe al de una población 
del principado de Mónaco, celebre por su casa de juego. 
Resulta que uno de los aparatos mecánicos más sencillos 
que permite obtener variables aleatorias es ... la ruleta. 
Acerca de esto hablaremos en el $ 3. En cambio, conviene, 
por lo visto, responder aquí a una pregunta frecuente: 
«¿Ayuda o no el método de Montecarlo a ganar en el juego 
de la ruleta?» No, no ayuda. E incluso no tiene nada que 
ver con este juego. 

1.2. Ejemplo. Para que el lector se haga una idea más 
clara de las cuestiones que vamos a tratar, consideremos 
un ejemplo muy sencillo. Supongamos que debemos calcular 

1) Metropolis N., Ulam S., The Monte Carlo me- 
thod, J. Amer. statistical assoc., 1949, 44, Ni 247, 335341. 

?) Se trata de los artículos de V. V. Chavchanidze, Yu, A. Sbréi- 

der y V. S. Vladímirov. 
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el área de una figura plana S. Ésta puede ser una figura 
arbitraria con una frontera curvilínea, definida gráfica o 
analíticamente y compuesta de uno o varios pedazos. Supon- 
gamos que se trata de la figura representada en la fig. 1 


viga 


y supongamos que toda la figura está comprendida dentro 
de un cuadrado de dimensión 4. 

Tomemos en el cuadrado un número N de puntos alea- 
torios. Sea N' el número de puntos que aparecen dentro de $. 
Es obvio, por razones geométricas, que el área de $ es igual 
aproximadamente a la razón N'/N. Cuanto mayor sea N 
tanto mayor será la exactitud de esta estimación. 

En el ejemplo representado en la fig. 1 se ha escogido un 
total de N = 40 puntos. De éstos, N” = 12 puntos aparecen 
dentro de S. Tenemos N'/N = 12/40 = 0,30 mientras que el 
valor exacto del área de $ es 0,35 1). 


1) En la práctica el método de Montecarlo no se 
aplica al cálculo de área de figuras planas ya que existen para ello 
otros métodos que, a pesar de ser más complejos, garantizan una 
exactitud mucho mayor. (Sigue en la pág. 11). 


1 


1.3. Dos peculiaridades del método de Montecarlo. 
La primera consiste en que su algoritmo tiene una estruc- 
tura muy sencilla. Como regla, se elabora primero un pro- 
grama para la realización de una prueba aleatoria (en el 
ejemplo del punto 1.2 hay que escoger un punto aleatorio 
en el cuadrado y comprobar si pertenece o no a 5). Después, 
esta prueba se repite V veces de modo que cada experimento 
sea independiente de los restantes y se toma la media de 
los resultados de todos los experimentos. 

Por esto el método de Montecarlo se denomina a veces 
método de pruebas estadísticas. 

La segunda peculiaridad consiste en que el error es, como 


regla, proporcional a la magnitud YD/N, donde D es una 
constante y N es el número de pruebas. Esta fórmula per- 
mite ver que para disminuir el error en 10 veces (en otras 
palabras, para obtener en el resultado otra cifra decimal 
exacta) es preciso aumentar N fo sea, el volumen del trabajo) 
en 100 veces. 

Se pone de manifiesto la imposibilidad de alcanzar en este 
camino una elevada exactitud. Por eso suele decirse que el 
método de Montecarlo resulta especialmente eficaz en la 
solución de problemas en les cuales se necesita conocer el 
resultado con poca exactitud (del 5 al 10 por ciento). 

Sin embargo, un mismo problema puede ser resuelto 
aplicando distintas variantes del método de Montecarlo *) 
a las que corresponden diferentes valores de D. En numero- 
sos problemas se logra elevar considerablemente la exacti- 
tud escogiendo un procedimiento de cálculo al que le co- 
rresponde un valor mucho menor de D. 


1.4. Problemas que permite resolver el método de Mon- 
tecarlo. En primer lugar, el método de Montecarlo permite 
simular cualquier proceso cuya marcha depende de factores 


Sin embargo, el método de Montecarlo expuesto en nuestro ejem- 
plo permite con la misma facilidad calcular el «volumen multidimen- 
sional» de un cuerpo en un espacio multidimensional. Es más, sucedo 
a menudo que el método de Montecarlo es en este caso el único método 
numérico que permite resolver el problema. 

2) En "la literatura especial con frecuencia se Suelo hablar aho- 
ra de los métodos de Montecarlo (en plural) para subrayar con 
ello que un mismo problema se puede resolver mediante la simulación 
de distintas variables aleatorias. 


» 
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aleatorios. En segundo lugar, en muchos problemas mate- 
máticos, que no tienen la menor relación con cuestiones alea- 
torias, se puede inventar un modelo probabilístico artificial 
(e incluso más de un modelo) que permite resolver estos 
problemas. En realidad, esto es lo que hemos hecho en el 
ejemplo del punto 1.2. 

Por consiguiente, se puede hablar del método de Mon- 
tecarlo como de un método universal en la solución de proble- 
mas matemáticos. 

Es de especial interés el hecho de que en algunos casos 
conviene abandonar la simulación del proceso aleatorio real 


y concentrarse en el análisis del modelo artificial. Una situa- 
ción de este tipo es considerada en el $ 7. 

1.5. Algo más sobre el ejemplo. Volvamos al ejemplo del 
punto 1.2. Para poder realizar el cálculo debemos escoger en 
el cuadrado puntos aleatorios, ¿Cómo realizar esto? 

Imaginemos el siguiente experimento. La fig. 1 (aumen- 
tada proporcionalmente), con la figura S y el cuadrado, se 
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coloca sobre un muro y sirve de blanco para un tirador que, 
situado a cierta distancia, dispara N veces apuntando al 
centro del cuadrado. Claro está que no todas las balas darán 
exactamente en el centro: perforarán el blanco en Y puntos 
aleatorios 1). ¿Es posible estimar el área de S a partir de 
estos puntos? 

En la fig. 2 viene representado el resultado de un experi- 
mento semejante. En este experimento se tiene N =40, 
N' =24 y N'/N =0,60 lo que representa casi el doble del 
valor exacto (0,35) del área. Por supuesto, tratándose de un 
tirador de gran puntería los resultados del experimento serán 
muy malos ya que casi todas las balas pegarán cerca del 
centro y perforarán $?). 

Es fácil persuadirse de que el método que hemos aplicado 
para calcular el área tendrá validez sólo si los puntos alea- 
torios, además de ser aleatorios, estarán «distribuidos uni- 
formemente» sobre el cuadrado. Para dar un sentido exacto 
a estas palabras es preciso estudiar la definición de variable 
aleatoria y algunas de sus propiedades. Estos conocimientos 
se exponen, sin demostración, en el $ 2. El lector que haya 
loído cualquier curso de la Teoría de las probabilidades puede 
omitir por completo el $ 2 a excepción del punto 2.5. 


1) Aceptamos que el tirador está lejos de ser el 
campeón del mundo y se encuentra a una distancia suficientemente 
grande del blanco. 

2) En el punto 4.5 explicaremos cómo han sido escogidos los 
puntos aleatorios en las figs. 1 y 2. 


CAPÍTULO 1 


SIMULACIÓN 
DE VARIABLES ALEATORIAS 


ye 
VARIABLES ALEATORIAS 


Suponemos que el lector está más o menos familiarizado 
con el concepto de probabilidad y pasamos directamente al 
concepto de variable aleatoria. 

Las palabras «variable aleatoria» suelen emplearse para 
subrayar que se ignora el valor concreto que tomará esta 
variable. Sucede también que estas palabras encierran inclu- 
so un desconocimiento de cuál es esta variable. 

En cambio, el matemático emplea estas mismas pala- 
bras «variable aleatoria» dándoles un contenido positivo 
bien determinado. Efectivamente, dice él, no conocemos 
el valor que tomará esta variable en un caso concreto dado 
pero sabemos qué valores puede tomar y también las proba- 
bilidades de unos u otros valores. A base de estos datos no 
podemos predecir con exactitud el resultado de una prueba 
relacionada con esta variable aleatoria pero sí podemos 
prever con gran seguridad los resultados de un gran número 
de pruebas. Cuanto mayor sea el número de pruebas, mayor 
exactitud tendrán nuestras predicciones. 

Es decir, para definir una variable aleatoria es preciso 
señalar los valores que puede tomar y las probabilidades 
de estos valores. 

2.1. Variables aleatorias discretas. Se dice que la variable 
aleatoría E es discreta si puede tomar un conjunto discreto 
de valores %,, La +... Tn”). 

Una variable aleatoria discreta E se define mediante la 


tabla 
Ly Lo... Tn ) 
= , E 
é e Pz +++ Pn (1) 
donde z,, %a, - - -, z, son los valores posibles de la variable 
E Y Pur Pas +» -» Pn Son las probabilidades que les corres- 


ponden. Hablando con más precisión, la probabilidad de 


1) La Teoría de las probabilidades estudia también 
variables aleatorias discretas que pueden tomar un número infinito 
de valores. 
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que la variable aleatoria £ tome el valor zx, (designaremos 
esta probabilidad por P (E =x,)) es igual a py: 


P(E=x) =p 


La tabla (T) se denomina distribución de la variable alea- 
toria, 
Hablando en términos generales, los números Z,, Tp, . + + 
. «+ Th pueden ser cualesquiera. En cambio, las probabili- 
dades P,, Pa, - - »» Pn deben cumplir dos condiciones: 
a) todos los números p, deben ser positivos: 


pi > 0; (0) 
b) la suma de todos los p; debe ser igual a 1: 
PM+FPDHRFP=1. (2) 


La última condición significa que E debe necesariamente 
tomar en cada caso uno de los valores %,, Za, .- ., Lp: 

Se denomina esperanza matemática de la variable aleatoria 
E el número. 


Mi= 3 Z1Pr> (3) 


Para aclarar el significado real de esta magnitud, la 
representaremos en la forma 


m=2 «nl 2, Pa 


De aquí puede verse que ME es el valor medio de la varia- 
ble E con la particularidad de que los valores x=, de mayor 
probabilidad figuran en la suma con pesos mayores *). 


2) La media ponderada aparece con frecuencia en 
las más diversas ramas de la ciencia. Por ejemplo, en la Mecánica: 
si en los puntos Z1, Xq, » +» » Zn (del eje Oz) están concentradas las ma- 
SAS Mi, May » » =y May la abscisa del centro de gravedad de esto sistema 
se dotermína mediante la fórmula 


Za sn 
3=Y 20m/ Y, mu. 
ini ii 


Por supuesto, la suma de todas las masas puede ser en este caso distinta 
de la unidad. 
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Señalemos las propiedades principales de la esperanza 
matemática: si c es una variable no aleatoria, se tiene 


M(E+0)=ME+<, (4) 

M (c5) =cME; (5) 

si E y n son dos variables aleatorias cualesquiera, se tiene 

M (E +m) = ME + Mn. (6) 

Se denomina varianza de la variable aleatoria £ el nú- 
mero 

D3 = M (E — M8)! (7) 


Por consiguiente, la varianza DE es la esperanza matemáti- 
ca del cuadrado de la desviación de la variable aleatoria E 
de su valor medio ME. Es obvio que siempre se tiene DE > 0. 

La esperanza matemática y la varianza son las caracterís- 
ticas numéricas más importantes de la variable aleatoria E, 
¿Qué valor práctico tienen? 

Supongamos que hemos observado varias veces la va- 
riable E obteniendo los valores E,, Ez, . . ., Ey (cada uno 
de los cuales es igual a uno de los números 2,, Ly, - +.» Zn); 


entonces la media aritmética de estos valores será próxima 
a ME 


L+b++.. HE) ME (8) 


y la varianza DE caracterizará el esparcimiento de estos 
valores alrededor del valor medio ME. 
Empleando (4), (5) y (6) podemos dar otra forma a (7): 


DE =M1£? — 2M5-E + ME)! =M (E) — 
— 2ME-M5 + ME)”, 


DE =M (E) — (ME. (9) 
Resulta más fácil, como regla, calcular la varianza apli- 
cando la fórmula (9) en lugar de la fórmula (7). 
Señalemos las propiedades principales de la varianza: 
si c es una variable no aleatoria, se tiene 
D (E +) =D8, (10) 
D (có) = c*DE- (11) 


de donde 
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En la Teoría de las probabilidades desempeña un papel 
importante el concepto de independencia de variables alea- 
torias. En realidad se trata de un concepto bastante com- 
plejo; pero en los casos más sencillos la cuestión se aclara 
fácilmente *). 

Para dos variables aleatorias independientes E y m son 
válidas las relaciones 

M (En) = ME-Mn, (12) 
D (5 + 1) = D5 + Dn. (13) 

Ejemplo. Consideremos la variable aleatoria x con la distribu- 

ción siguiente 
1234556 
(8 vs 4/6 4/6 4/6 1/6)" 


Es obvio que podemos realizar esta variable considerando el número de 
puntos que aparece al arrojar un dado: todos los valores son igual- 
mente probables. Calculemos la esperanza matemática y la varianza 
de x. Según la fórmula (3), tenemos 


4 1 
Milpa co +63, 
Según la fórmula (9), tenemos 
Dx=M (x?) —(Mx)2= n+ co (8,5% =2,917. 


Ejemplo, Consideremos la variable aleatoria O con la distribución 
siguiente 
o=( 34 ) 
la a, 


Podemos realizar esta variable considerando el juego de cara y cruz en 
el quo a la cara lo corresponden 3 puntos y a la cruz 4 puntos. En este 
caso tenemos 
MO = 0,5-3 + 0,5:4 == 3,5; 
DO = 0,5 (3? + 42) — (3,5)* = 0,25, 
Como vemos MO = Mx mientras que DG < Dx, Lo último era de 


cp ya que los valores de € pueden diferir de 3,5 en -+ 0,5 todo lo 
más mientras que la dispersión de los valores de x puedo alcanzar 


hasta + 2,5. 


2.2. Variables aleatorias continuas. Supongamos que 
en el origen de coordenadas del plano se encuentra una can- 


1) Supongamos que además de la variable E estamos 
observaudo otra variable aleatoria y. Si la distribución de la variable E 
no se altera al conocerse el valor que ha tomado la variable y, es 
lógico decir que E no depende de y. 
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tidad determinada de radio. Todo átomo de radio emite 
al desintegrarse una partícula a. Su dirección puede ser 
caracterizada mediante el ángulo y (fig. 3). Tanto desde 


el punto de vista teórico como práctico, la emisión es po 
sible en cualquier dirección; por eso, esta variable alea- 
E puede tomar cualquier valor comprendido entre 0 
y 2. 

Diremos que una variable aleatoría E es continua si puede 
tomar cualquier valor comprendido en un intervalo (a, b). 

Toda variable aleatoria continua E queda definida si 
se da el intervalo (a, b) que contiene los valores posibles 
de esta variable y la función p (x) que lleva el nombre de 
densidad de probabilidad de la variable aleatoria E (o den- 
sidad de distribución de E). 

El significado real de la función p (x) es el siguiente: 
sea (a', b”) un intervalo cualquiera contenido en (a, b) (es 
decir, sea a <a' y b'< b); entonces la probabilidad de 
que É tome un valor perteneciente al intervalo (a”, b”) 
es igual a la integral 

pr 
P (a <E<0)=Í pla) dz. (14) 


ES 


20 


El área sombreada en la fig. 4 es igual al valor de la inte- 


gral (14). 
El conjunto de los valores de £ puede formar un intervalo 
cualquiera. Puede darse incluso el caso de que a = —oo 


y también el de que b= oo. En cambio, la densidad p (x) 
debe cumplir dos condiciones análogas a las condiciones (1) 
y (2) para las variables discretas: 

a) la densidad p (x) debe ser positiva: 


p (2) >0; (15) 
b) la integral de la densidad p (x) correspondiente a todo 
el intervalo (a, b) debe ser igual a 1: 


b 
$ pla) dr=1. (16) 


z 
Se denomina esperanza matemática de la variable aleato- 


ria continua el número 


b 
Mi= j zp (2) dz. an) 


El significado de esta característica es el mismo que en 
el caso de una variable aleatoria discreta. Efectivamente, 
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puesto que d ó 
MiE= $ xp(x) 22] $ p(x)dz, 


es fácil ver que se trata del valor medio de E pues cualquier 
valor de É es un número zx del intervalo (a, b) que aparece 
en la integral con el peso p (x) *). 

Todo lo expuesto en el punto 2.1, desde la fórmula (4) 
hasta la fórmula (13) inclusive, sigue siendo válido también 
para las variables aleatorias continuas; esto se refiere a la 
definición (7) de la varianza, a la fórmula (9) que se emplea 
para calcularla y a todas las propiedades de ME y DE. La 
repitición huelga. 

Mencionemos únicamente otra fórmula más que permite 
calcular la esperanza matemática de una función de E. 
Sea, como antes, p (z) la densidad de probabilidad de la 
variable aleatoria E. Tomemos una función continua cual- 
quiera f (x) y consideremos la variable aleatoria y =1 (E). 
Se puede demostrar que 


1 
MO=) 1()p (0 dz. (18) 


Cabe subrayar que, hablando en términos generales, 
M/ e 1 (ME). 

i la variable aleatoría y está definida en el intervalo 
(0, 1) y su densidad es p (x) =1 (fig. 5), se dice que está 
uniformemente distribuida en (0, 1). 

Efectivamente, cualquiera que sea el intervalo (a”, b') 
que se tome dentro de (0, 1), la probabilidad de que y tome 
un valor perteneciente al intervalo (a', b”) es igual a 


3) También en este caso podemos señalar una fór- 
mula análoga de la Mecánica: si p (z) es la densidad lineal de la barra 
«<< b, la abscisa del centro de gravedad de la barra so calcula 
mediante la fórmula 


d b 
z= $ ap (m1 02 $ p(=) dz. 
a 


a 
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o sea, es igual a la longitud de este intervalo. En particular, 
si dividimos (0, 1) en un número arbitrario de intervalos 


P(x) 


Fig. 5 


de igual longitud, la probabilidad de que y tome un valor 
perteneciente a cualquiera de estos intervalos será la misma. 
Es fácil calcular que 


My= 


1 
xp (2)dx= Ñ ado 
0 


4 


Dy= | 2%p (2) d2— (MY == =+3 + 


1 
0 
1 

La variable aleatoria y será empleada frecuentemente 
en lo sucesivo. 

2.3. Variables aleatorias normales. Se da el nombre de 
variable aleatoria normal (o gaussiana) a toda variable 
aleatoria E que está definida en todo el eje (—oo, 00) y que 
tiene la densidad 
(x-aj2 


ear, (19) 


L 

1) = == 
p(=) VEZ 
donde a y o >0 son unos parámetros numéricos Dl 


2) o es un número y no una variable aleatoria. Sin 


embargo, ompleamos esta lotra griega ya por tradición. 
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El parámetro a no influye en la forma de la curva p (2): 
su variación sólo conduce a un desplazamiento de la curva a 
lo largo del eje x. En cambio, al variar O sí se altera la forma 
de la curva. En efecto, es fácil ver que 

E: 
máx p (2) =p(0)= 5: 

Si disminuimos o, aumentará el máx p (x). Pero toda el 
área que encierra la curva p (x) es igual a 1 según (16). Por 


5=0,5 


-3 -2 Y o ] 2 3 


Fig. 6 


eso, la curva tendrá que estirarse hacia arriba en una vecin- 
dad de x =au y decrecer para todos los valores suficiente- 
mente grandes de z. En la fig. 6 hemos representado dos 
densidades normales correspondientes al caso a=0 y 
oc=1 y al caso a =0 y o =0,5. (Otra densidad normal 
será construida en la fig. 21 de la página 54). 

Se puede demostrar que 


Mí =a y DE =0?. 


Las variables aleatorias normales aparecen con gran 
frecuencia al analizar cuestiones de la más diversa índole. 
La razón de esto será explicada a continuación. Por ejemplo, 


2 


el error de medición $ es, como regla, una variable aleatoria 
normal. Si no existe error sistemático en las mediciones, se 
tiene a=M6=0. La magnitud o =YD8, llamada 
error medio cuadrático, caracteriza el error del método que 
se emplea para la medición. 

La regla de «las tres sigmas». Es fácil calcular que 


0430 
Pp (2) dx =0,997 
aldo 
cualesquiera que sean a y o en (19). Do (14) resulta que 
P(a—30<t<a+ 30) = 0,997. (20) 


La probabilidad 0,997 es tan próxima a 1 que la última 
fórmula se interpreta a veces así: al efectuar una prueba es 
prácticamente imposible obtener un valor de L que difiera de 
Mi en más de 30. 

2.4. Teorema central del límite en la Teoría de las proba- 
bilidades. Laplace fue el primero en enunciar este importante 
teorema. Varios matemáticos destacados, entre los cuales 
figuran P. L. Chébishev, A. A. Márkov y A. M. Liapunoy, 
estudiaron sus generalizaciones. Su demostración es bastante 
compleja. 

Consideremos N variables aleatorias E,, En, ..-., En 
independientes e idénticas de modo que las distribuciones de 
probabilidad de estas variables coinciden. Por consiguien- 
te, también coinciden sus esperanzas matemáticas y varian- 
zas. 

Pongamos 


ME, = Mé, =... =Mióy =m, 
DE =D3=... = Diy = DB 
Designemos por py la suma de todas estas variables: 
y = E 4 Es + Eno 
De las fórmulas (6) y (13) se deduce que 
Moeyn =M(E + E. +...+ En) = Nm, 
Dey =D(E-+E2+ +... + En) = NL 


Consideremos ahora la variable aleatoría normal £y 
con estos mismos parámetros: a = Nm y 07 = NO, 
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El teorema central del límite afirma que cualquiera que 
sea el intervalo la”, b”) se tiene para valores grandes de N 


e 
Pl <ox <0) 2 Í pr, (2) dz. 
2 


El significado real de este teorema es obvio: la suma py 
de una gran cantidad de variables aleatorias idénticas es 
aproximadamente normal (pp, (2) == Pr, (2))- 

En realidad, este teorema es válido en situaciones mucho 
más generales: los sumandos E, Ez, . . «, Ey pueden no ser 
idénticos e independientes; lo que importa es sólo que cada 
sumando por separado no influya considerablemente en la 
suma. 

Es este teorema el que explica porque las variables alea- 
torias normales aparecen con tanta frecuencia en la práctica. 
En efecto, siempre que tropezamos con el efecto sumario de 
una gran cantidad de factores aleatorios insignificantes, la 
variable aleatoria resultante es normal, 

Por ejemplo, la desviación del proyectil respecto al obje- 
tivo es casi siempre una variable aleatoria normal por cuanto 
depende de las condiciones meteorológicas en las diferentes 
partes de la trayectoria y de otros muchos factores. 

2.5. Esquema general de aplicación del método de Monte- 
carlo. Supongamos que debemos calcular una magnitud m 
que desconocemos. Tratemos de idear una variable aleatoria £ 
tal que MÉ = m. Sea además DÉ = b?, 

Consideremos N variables aleatorias independientes É,, 
En - + -, Ey con la misma distribución que tiene E. Si 
NÑ es suficientemente grande, la distribución de la suma 
er =E+F EA. -. + Ey será, según el teorema del 
punto 2.4, una distribución aproximadamente normal con 
los parámetros a = Nm y 0? = Nb?. De (20) se deduce que 


P (Vm—3VN< py < Nm + 30Y N) = 0,997. 
Si dividimos entre N la desigualdad que figura entre las 
llaves, obtendremos una desigualdad equivalente; la proba- 
bilidad continuará siendo la misma: 


P (m7 <<m+ 77) 0,997. 


3-1216 
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Podemos representar este último resultado en la forma Si- 
guiente 


N 
P (| 35—m|<- 77) = 0,997. (21) 
jui 


Ésta es una relación de suma importancia para el método 
de Montecarlo. Ofrece un método para calcular m y a la vez 
permite estimar el error. 

En efecto, determinemos N valores de la variable aleato- 
ria E 2). Según la fórmula (21), la media aritmética de estos 
valores será igual aproximadamente a m. Con una elevada 
probabilidad podemos afirmar que el error de esta aproxima- 
ción no pasa de 3b/V N. Es obvio que el error tiende a cero 
cuando N crece. 


$3. 
OBTENCIÓN 
DE VARIABLES ALEATORIAS 
EN LAS MCE 


El propio planteamiento de la cuestión («obtención de 
variables aleatorias en las MC E») provoca cierta perplejidad: 
¿de qué casualidad se puede hablar si todo cuanto realiza la 
máquina debe ser programado de antemano? 

En esta cuestión hay efectivamente ciertas sutilezas; 
pero conciernen más bien a la filosofía y, por eso, las deja- 
remos a un lado. 

Sin embargo, cabe señalar que las variables aleatorias 
consideradas en el $ 2 representan conceptos matemáticos 
ideales. Sólo la práctica dilucida si pueden ser empleados 
para describir determinados fenómenos de la naturaleza. 
Esta descripción es siempre aproximada. Es más, una va- 
riable aleatoria que describe satisfactoriamente una magnitud 
física en unos fenómenos, puede resultar una característica 


1) Resulta lo mismo determinar un valor de cada 
una de las variables Es, Ea, . » ., Ey o hallar N valores de la variable E 
ya quo todas estas variables aleatorias son idénticas (tienen la misma 
distribución). 
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desafortunada de esta misma magnitud al tratarse de otros 
fenómenos. 

De la misma forma, un camino que se asemeja a una 
recta (a una recta matemática ideal «sin anchura») en el 
mapa del país, se convierte en una franja sinuosa en el 
plano a gran escala de una localidad... 

Suelen distinguirse tres procedimientos que se aplican 
para obtener variables aleatorias: tablas de números alea- 
torios, generadores de números aleatorios y el método de 
números seudoaleatorios. 

3.1. Tablas de números aleatorios. Realicemos el experí- 
mento siguiente. Tomemos diez tarjetas iguales escribiendo 
en ellas las cifras 0, 1, 2, . . ., 9, metamos todas las tarjetas 
en un sombrero y, después de mezclarlas bien, comencemos a 
extraeslas, con la particularidad de que cada vez sacaremos 
una tarjeta solamento, la volveremos a meler on el sombrero 
y mezclaremos de nuevo todas las tarjotas. Las cifras así 
obtenidas pueden ser recogidas en una tabla somejante a la 
tabla A de la página 78 al final del libro (por razones de 
comodidad las cifras en la tabla A vienen agrupadas de cinco 
en cinco). 

Esta tabla lleva el nombre de tabla de números aleatorios 
aunque sería más justo denomiuarla tabla de cifras aleato- 
rias. Dicha tabla puede ser introducida en la memoria de la 
MCE. Cada vez que, durante los cálculos, necesitemos un 
valor de la variable aleatoria con la distribución 


a 
laa 0,1 0,1 AJA (22) 


podemos tomar una tras otra las cifras de esta tabla. 

La mayor tabla de números aleatorios publicada hasta el 
momento contiene 1 000 000 de cifras!). Por supuesto, fue 
obtenida empleando, en lugar del sombrero, una técnica más 
moderna: una ruleta electrónica. Un esquema elemental de 
esta ruleta aparece en la fig. 7 (el disco giratorio es para- 
do en seco y se toma la cifra que señala la flecha). 

Es preciso subrayar que no es tan fácil como parece obte- 
ner una tabla buena de números aleatorios. Todo aparato 


y 1) RAND Corporation, A atillion random digits 
with 100 000 normal deviates, The Free Press, 1955. 


y 
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físico real genera variables aleatorias cuyas distribuciones 
difieren ligeramente de la distribución ideal (22). Además, 
durante el experimento pueden surgir errores (por ejemplo, 
una de las tarjetas puede adherirse por cierto tiempo al 
forro). Por eso, las tablas obtenidas suelen comprobarse 


Fig. 7 


minuciosamente a base de tests estadísticos especiales para 
ver si no hay contradicción entre unas u otras propiedados 
de grupos de números y la hipótesis de que estos números 
representan valores de la variable aleatoria (22). 

Veamos uno de los tests más sencillos. Consideremos una 
tabla que contiene N cifras. Supongamos que el número de 
ceros en esta tabla es v,, el número de unos es v,, el número 
de doses es v,, etc. Calculemos la suma 


9 
2 (v¿—0,1)2 


La Teoría de las probabilidades permite predecir los 
límites entre los cuales puede estar comprendida esta suma 
que no puede ser excesivamente grande (por cuanto la espe- 
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ranza matemática de cada v, es igual a 0,1 N) ni excesiva- 
mente pequeña (ya que esto equivaldría a una «regularidad 
excesiva» en la distribución de los valores). 

Las tablas de números aleatorios se emplean sólo cuando 
los cálculos correspondientes al método de Montecarlo se 
realizan a mano. Esto se debe a que todas las MCE poseen una 
memoria interna relativamente reducida en la que no se 
puede introducir una tabla grande. Por otra parte, tampoco 
conviene introducir la tabla en la memoria externa y recu- 
rrir constantementea ella ya que esto afectaría considerable- 
mente la velocidad del cálculo. 

Sin embargo, no está excluido que, con el tiempo, la 
capacidad de memoria de las MCE aumente sensiblemente 
y entonces las tablas de números aleatorios resulten de gran 
utilidad. 

3.2. Generadores de números aleatorios. La ruleta mencio- 
nada en el punto 3.1 podría, al parecer, combinarse con la 
computadora y emplearse para generar, según sea necesario, 
números aleatorios. Pero todo dispositivo mecánico resulta 
demasiado lento comparado con la MCE. Por eso, para 
generar variables aleatorias suelen emplearse los ruidos 
en las lámparas electrónicas: si durante un intervalo fijo At 
de tiempo el nivel del ruido sobrepasa el umbral escogido 
un número par de veces, se inscribe el cero; si esto ocurre 
un número impar de veces, se inscribe el uno?). 

A primera vista, éste es un procedimiento muy cómodo. 
Imaginemos m generadores de este tipo que funcionan sin ce- 
sar paralelamente y envían ceros y unos aleatorios a todos 
los órdenes binarios de una célula especial. Cada ciclo co- 
rresponderá a un número de m órdenes. Durante el cálculo 
podemos dirigirnos en cualquier momento a esta célula y 
tomar de ella un valor de la variable aleatoria y uniforme- 
mente distribuida en el intervalo (0, 1). Este valor, claro 
está, será aproximado y quedará representado por una frac- 
ción binaria de m órdenes 0, %y% ... Gm donde cada 
Q. imita la variable aleatoria con la distribución 


(o. 1 
Y 


») Existen procedimientos más perfectos, 
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Sin embargo, también este método tiene sus defectos. 
En primer lugar, es difícil comprobar la «calidad» de los 
números generados. Se hace necesario un control periódico 
ya que debido a distintos desarreglos puede surgir la así 
llamada «deriva en la distribución» (cuando los ceros y los 
unos en uno de los órdenes dejan de aparecer con la misma 
frecuencia). En segundo lugar, todos los cálculos en las 
MCE suelen realizarse dos veces para evitar intermitencias 
casuales, Pero, sin recurrir durante los cálculos a la memo- 
ria, es imposible reproducir los mismos números aleatorios 
y si recurrimos a la memoria, volvemos al caso de las tablas. 

Generadores de este tipo resultarán útiles, sin duda, 
cuando comenzarán a construirse MCE especializadas, adap- 
tadas a la solución de problemas por el método de Monte- 
carlo. En las MCE universales, que raramente se aplican 
para los cálculos a base de números aleatorios, no es econó- 
mico mantener y explotar un dispositivo especial. Resulta 
más conveniente emplear los así llamados números seudoalea- 
torios. 

3.3. Números seudoaleatorios. Como quiera que la «cali- 
dad» de los números aleatorios empleados se controla median- 
te tests especiales, podemos no interesarnos por el procedi- 
miento aplicado para obtenerlos exigiendo solamente que 
verifiquen el sistema de tests adoptado. Incluso podemos 
intentar calcularlos a partir de una fórmula dada. Por su- 
puesto, debe ser una fórmula muy ingeniosa. 

Los números que.se obtienen a partir de una fórmula y que 
imitan los valores de la variable aleatoria y se denominan 
números seudoaleatorios. La palabra «imitan» significa que 
estos números verifican una serie de tests como si fuesen 
valores de esta variable aleatoria. 

El primer algoritmo destinado a la construcción de núme- 
ros seudoaleatorios fue propuesto por J. von Neumann y se 
conoce como el método de centros de cuadrados. Para explicarlo 
recurriremos a un ejemplo. 

Consideremos el número y. = 0,9876 formado por cuatro 
cifras. Al elevarlo al cuadrado obtendremos el número 
y! = 0,97535376 de ocho cifras. Tomemos las cuatro cifras 
que aparecen en el centro de este número y pongamos y, = 
= 0,5353. 

Elevando ahora y, al cuadrado (yí = 0,28654609), escoja- 
mos de nuevo las cuatro cifras del centro y pongamos Ya = 
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= 0,6546. Procediendo de la misma forma, tendremos y = 
= 0,42850116, y, = 0,8501; y? = 0,72267001, y, = 0,2670; 
y = 0,07128900, ys = 0,1289, etc.!). 

Pero este método no puede considerarse satisfactorio pues 
da más números pequeños de lo que se necesita. Por eso, 
fueron desarrollados por distintos científicos otros algoritmos. 
Algunos se basan en las particularidados de las MCE con- 
cretas. A título de ejemplo, consideremos uno de estos algo- 
ritmos empleado en la MCE Strelá. 


Ejemplo, 2). Strelá es una MCE de tres direcciones con coma flolan- 
te. La célula en la que se inscribe el número z consta de 43 órdenes bina- 
rios (fig. 8). La computadora opera con números binarios en forma 
normalizada: == + 4-2*P, dundo p es el orden y q la mantisa del 


t E 
Signo de montisa Signo del orden 


Fig. 8 


número. En el orden j-ésimo de la célula puede figurar un cero o un 
uno; representemos por ej esta variable. Entonces 
e, y 95 
ir ad El 


P=e3125+e324+... 421920. 


La condición de normalización 0,5 < y < 1 significa que e, es siem- 
pre igual a 1. El signo < + > se representa por el cero y el signo <-> 
por el uno. o 
. El número yy4, se obtiene del número yy realizando tres opera- 
ciones: 

4) ya se A ica por una constanto grande que, como regla, 
so toma igual a 4 


2) Podemos resumir este algoritmo en la fórmyla 
Yaga =P (1), dondo F representa el comjuuto de operaciones que 
debemos realizar con el número y, para obtener ya41 El número yo eS 
un número dado. 

2) Véase H. M. Co6oab, Nconnocayualimo wueno MIE MAnEN 
Cmpesa, Teopma seponuocri 1 ee npuaenemns, 3, Ni 2 (1958), 205. 
244 (1. M- Sóbol, Números soudoaleatorios para la computadora Strelá; 
en la revista «Teoría de las probabilidados y sus aplicaciones», 3, 
n* 2 (1958), pág. de 205 a 241). 
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2) la imagen del producto 101” y; se desplaza en siete órdenes hacia 
la izquierda (de modo que los siete primeros órdenes del producto desa- 
parecen mientras que en los órdenes del 36 al 42 aparecen ceros); 

3) se toma el valor absoluto del número obtenido (después de nor- 
marlo) que será precisamente Yp+r+ 

Si arrancamos de yo == 1, este proceso permite obtener más de 
80 000 múmeros diferentes ya; después en la sucesión surgo un período 


y los números comienzan a repetirse, Diferentes controles de los 
50 000 primeros números han dado resultados plonamento satisfacto- 
rios, Estos números se hau empleado más de una vez para la solución 
de diversos problemas. 


Las ventajas del método de números seudoaleatorios 
son evidentes. En primer lugar, para obtener un número 
basta realizar unas cuantas operaciones sencillas y gracias 
a esto la velocidad de generación de números aleatorios es 
del mismo orden que la velocidad de funcionamiento de la 
MCE. En segundo lugar, el programa ocupa sólo unas cuan- 
tas células de la memoria, En tercer lugar, cualquiera de 
los números y» puede ser reproducido fácilmente. Por últi- 
mo, basta comprobar una vez solamente la «calidad» de esta 
sucesión para poder aplicarla después reiteradamente y con 
seguridad en la solución de problemas semejantes. 

El método tiene sólo un defecto: la «reserva» de los nú- 
meros seudoaleatorios es limitada. Pero existen procedi- 
mientos que permiten obtener una cantidad mucho mayor 
de números. En particular, se pueden variar los números 
inicialos ys. Actualmente los cálculos en los que se emplea 
el método de Montecarlo se realizan, en su aplastante mayo- 
ría, aprovechando los números sendoaleatorios. (Véanse in- 
formaciones complementarias en le $ 10). 


54 
TRANSFORMACIONES 
DE VARIABLES ALEATORIAS 


En la solución de diferentes problemas tropezamos con la 
necesidad de simular distintas variables aleatorias. En las 
etapas iniciales del desarrollo del método de Montecarlo los 
investigadores intentaban obtener cada variable aleatoria 
construyendo una ruleta apropiada. Por ejemplo, para hallar 
los valores de la variable aleatoria con la distribución 
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siguiente 


HQ E En . 
los 0,25 0,125 El (22%) 


se puede emplear la ruleta representada en la fig. 9 que fun- 
ciona igual que la ruleta de la fig. 7, pero no está dividida 
uniformemente sino en partes proporcionales a las py. 


Fig. 9 


Pero esto no hace falta: resulta que los valores de cual- 
quier variable aleatoria se pueden obtener efectuando trans- 
formaciones de una sola variable aleatoria (es decir, de una 
variable estandard»). Suele emplearse con este fin la variable 
aleatoria y uniformemente distribuida en el intervalo (0, 1). 
Conocemos ya cómo se obtienen los valores de y. 

Denominaremos sorteo de la variable aleatoria E el pro- 
ceso de determinación de un valor de esta variable aleatoria E 
mediante la transformación de uno o varios valores de y. 

4.1. Sorteo de una variable aleatoria discreta. Suponga- 
mos que es preciso obtener valores de la variable aleatoria E 


con la distribución 


E e 
s- (7. ¿Pa 


Consideremos el intervalo 0 < y < 1 dividiéndolo en nr 
intervalos de longitudes p,, Pa, » + -; Pn- Es evidente que 
las coordenadas de los puntos de división serán y <= p,, 
Y=P HP» Y=PM+P2 HP =P + P4 
“+... Pa-1> Indiquemos con las cifras 1, 2, ..., n los 


——————_ z———_—_—_—e _  —_u—— 
0 PL Pipa TN 1-Pn y 
Fig. 10 


intervalos obtenidos (fig. 10). Con esto concluyen los prepa- 
rativos del sorteo de E. Cada vez que tengamos que «reali- 
zar un experimento» y sortear el valor de E, tomaremos un 
valor de y y construiremos el punto y = y. Si este punto 
aparece en el intervalo correspondiente al número ¿, acep- 
taremos que E = z, (en este experimento). 

Es muy fácil argumentar la legitimidad de este procedi- 
miento. En efecto, puesto que la variable aleatoria y está 
uniformemente distribuida en (0, 1), la probabilidad de que 
y pertenezca a uno de los intervalos es igual a la longitud 
del mismo. Por consiguiente, 


P(0<7<p)=Pps 
Pp <Y< pi +2) = Pa, 


P (DH Pad o. Pas < Y <1) = Pro 
Según nuestro método, tomamos É = x, cuando 
PEPE PA PSP Pa pe 


y la probabilidad de ello es igual a py. 

Por supuesto, al emplear la MCE podemos pasarnos sin 
la fig. 10. Supongamos que tanto los números 2, Za, - + ., Ln 
como las probabilidades P,, P, + Pa, Py + P2 + P3 + .., 1 


35 


han sido introducidos uno tras otro en la memoria. En la 
fig. 11 de la pág. siguiente representamos el esquema 
sinóptico del subprograma de sorteo de E. 


Ejemplo. Soriéense 10 valores de la variable aleatoria con la dis- 
tribución 


Us (esa dea) . 


Escojamos los valoros de y tomando pares de cifras de la tabla A 
(pág. 78) multiplicados por 0,01 1). Es decir, y — 0,86; 0,54; 0,59; 
0,07; 0,95; 0,68: 0,15; 0,56; 0, ,34. 

Según nuestro esquema, los valores de y menores que 0,58 corres 
ponden a U = 3 mientras que los valores y > 0,58 corresponden a 0 = 
= 4; de modo que obtenemos los valores siguientes: Y — 4 4; 
4; 3,3; 4,3, 

Cabe subrayar que el orden de numeración de los valores 
Zi Lg - - -, % en la distribución de E puede ser arbitrario 
pero debe fijarse antes de comenzar el sorteo. 

4.2. Sorteo de una variable aleatoria continua, Supone 
mos ahora que debemos obtener valores de una variable 
aleatoria ¿ distribuida con la densidad p (x) en el interva- 
lo (a, b). 

Demostremos que los valores de E so pneden determinar 
de la ecuación 


¡ p(a)dr=y, (23) 


a 


o sea, que escogido ya el valor do y es preciso revolver la 
ecuación (23) para hallar el valor de E. 
Para demostrar esta afirmación consideremos la función 


ya $ P() dz. 


1) Puesto que on este ejemplo los valores de p, tienen 
dos cifras decimales, basta tomar los valores de y con dos cifras 
decimales. Aceptado este grado de aproximación, so nos puede dar el 
caso de y = 0,58 que debe ser agregado al caso de y > 0,58 (ya quo 
el valor y=0,00 es posible mientras que el valor y ,00 es 
imposible). Si empleamos para y números con varias cifras decimales, el 
caso y=p, se hace poco probable y puede ser agregado a cualquiera 
de las desigualdades. 
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Determinamos y 


] 
l 
== 
(A 


Aumentamos en 1 las direc= 
ciones de las células de las 
que se toman Py y 21 


Reconstruimos las 
direcciones 


de pi y za 


Fig. 11 
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De las propiedades generales (15) y (16) do la densidad resulta 
que 


y (a) =0, y (0) =1 
y que la derivada 
y' (x) =p (2) >0. 


Por lo tanto, la función y (x) crece monótonamente del O 
al 1 (fig. 12). Debido a esto, toda recta y = y, donde 0 < 
<wy< lt, corta el gráfico de y = y (1) en un punto único 


Fig. 12 


cuya abscisa se toma precisamente igual a E. Es decir, la 
ecuación (23) tiene siempre solución y esta solución es única. 

Tomemos ahora un intervalo arbitrario (a”, b') conte- 
nido en (a, b). A los puntos 


a<z<b 
de este intervalo les corresponden las ordenadas de la curva 
y =y(z) que satisfacen la desigualdad 
y (a) <y< y (0). 


Por eso si E pertenece al intervalo a” <z< bh”, resulta 
que y pertenece al intervalo y (a) < y < y (b') y viceversa 


(fig. 13). Luego 
Pl <E<D) =P (y (a) <y< y (6). 


Puesto que y está uniformemente distribuida en (0, 1), 
tenemos 


A 
Pa <r<r 0) =4(0)—y (e) = | p (a) dz. 
Es decir, 


9 
Pla <E<0)=f p (2) dz; 


pero esto significa precisamente que la variable aleatoria 


o o ay b z 
Fig. 13 
— A A  —_ 


E, raíz de la ecuación (23), tiene la densidad de probabilidad 
Pp (2). 


Ejemplo, Se dice que la variable aleatoria n está uniformemente dis 
eres en el intervalo (a, b) si su densidad es constante en este inter- 
valo: 


1 
p()=5= Para a<=<b, 
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Para sortear los valores de y formamos la ecuación (23): 


La integral se calcula fácilmente y resulta 


al 
EAN 


do donde obtenemos una expresión explícita para n: 
n=2e+ y(b—a). (24) 


Otros ejemplos de aplicación de la fórmula (23) pueden 
verse en los puntos 5.2 y 8.3. 

4.3. Método de Neumann para el sorteo de una variable 
aleatoria continua. Puede resultar que sea muy difícil resol- 
ver la ecuación (23) respecto a E; por ejemplo, si la integral 
de p (x) no se expresa en funciones elementales o si la den- 
sidad p (zx) viene dada gráficamente. 

Supongamos que la variable aleatoria E está definida 
e O finito (a, b) y que su densidad está acotada 
(fig. , 


Pp (2) < Mo. 


Los valores de ¿ se pueden sortear de la siguiente manera 

1) se escogen dos valores y” y y'* de la variable aleatoria y 
y se considera el punto aleatorio T' (n', 1'*) con las coor- 
denadas 


v=a+ Y (ba y0" =y9Mu 


2) si el punto T' aparece por debajo de la curva y = 
= p (2), tomamos E = n'; si el punto T' aparece por encima 
de la curva y =p (x), retiramos de la consideración el 
par (y, y'”) y escogemos un par muevo de valores (y”, y””). 

La argumentación de este método viene en el punto 9.1. 

4.4, Sobre el sorteo de variables normales. Existe una 
gran variedad de procedimientos que se aplican para sortear 
distintas variables aleatorias. No nos detendremos en ello. 
Suelen emplearse sólo si los métodos de los puntos 4.2 y 4.3 
resultan poco efectivos. 
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En particular, esto ocurre en el caso de la variable alea- 
toria normal £ ya que la ecuación 


j Faz =7 


e e 
VE 


no admite solución explícita y el intervalo de los valores posi- 
bles de £ es infinito. 


Fig. 14 


La tabla B que viene al final del libro (en la pág. 78) 
contiene valores (ya sorteados) de la variable aleatoria nor- 
mal Z con la esperanza matemática Mí =0 y la varianza 
Dí = 1. Es fácil demostrar!) que la variable aleatoria 


E =a+of (25) 
también será normal; además, de (10) y (11) se deduce que 
MÉ' =a y DE' =0?. 


1) La demostración viene en ol punto 9.2 
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Por consiguiente, la fórmula (25) permite, empleando la 
tabla B, sortear valores de cualesquiera variables normales. 

4.5. Otra vez sobre el ejemplo del punto 1.2. Ahora pode- 
mos explicar cómo se han escogido los puntos aleatorios en 
las figuras 1 y 2. En la fig. 1 hemos construido los puntos 
con las coordenadas 

r=yey=py”. 

Los valores y” y y'* se han calculado a base de los grupos 
de cinco cifras de la tabla A: x, = 0,86515, y, = 0,90795; 
z, =0,66155, y, = 0,66434, etc. 

Se puede demostrar*) que, debido a que las abscisas y las 
ordenadas de nuestros puntos son independientes, la probabi- 
lidad de que uno de estos puntos aparezca en una región 
cualquiera interior del cuadrado es igual al área de esta 
región. En otras palabras, los puntos están uniformemente 
distribuidos en el cuadrado. 

4 En la fig. 2 hemos construido los puntos con las coordena- 
as 
zx =0,5 + 0,22 e y =0,5 + 0,28'* 


con la particularidad de que los valores Z' y ['” se toman 
de la tabla B uno tras otro: 


x, = 0,5 + 0,2-0,2005, y, = 0,5 + 0,2 -1,1922; 
2, =0,5 + 0,2(—0,0077), ... 
Un punto no se ha considerado por quedar fuera del cuadra- 
do. 
De la fórmula (25) se deduce que las abscisas y las ordena- 
das de estos puntos representan variables aleatorias normales 
con la esperanza matemática a = 0,5 y la varianza 0* = 0,04 


*) La demostración viene en el punto 9.3. 


CAPITULO 2 
EJEMPLOS 
DE APLICACIÓN 
DEL MÉTODO 
DE MONTECARLO 


$5. 
ANALISIS 
DE UN SISTEMA DE SERVICIOS 


5.1. Descripción del problema. Consideremos uno de 
los sistemas más sencillos de servicios. Consta de n líneas 
(canales o centros de servicio), cada una de las cuales se 
puede emplear para «atender a los consumidores». El sistema 
recibe pedidos, con la particularidad de que los momentos 
de recepción son aleatorios. Todo pedido va a la línea núme- 
ro 1. Sien el momento de recepción del k-ésimo pedido (de- 
signemos este momento por 7,) la línea está libre, comienza 
a cumplir el pedido durante un tiempo too (toc €s el tiempo de 
ocupación de la línea). Si en el momento 7, la línea número 1 
está ocupada, el pedido se transmite instantáneamente a la 
línea número 2, etc. Por último, si en el momento 7, resul- 
tan ocupadas todas las líneas, el sistema rechaza el pedido. 

Es preciso determinar el promedio de pedidos que atenderá 
y rechazará el sistema durante un tiempo 7. 

Está claro que semejantes problemas pueden surgir al 
analizar el funcionamiento de cualquier empresa y no sólo 
de las empresas de servicios. En algunos casos muy espe- 
ciales se logra la solución analítica. Pero, en casos complejos 
(de los que trataremos más tarde), el método de Montecarlo 
constituye el único método de cálculo. 

5.2. Flujo elemental de pedidos. El primer problema que 
se plantea al analizar sistemas de este tipo es: ¿qué repre- 
senta el flujo de pedidos? Este problema se resuelve de mane- 
ra experimental a base de un estudio suficientemente pro- 
longado de los pedidos. El análisis de los flujos de pedidos, 
correspondientes a distintas condiciones, ha permitido desta- 
car algunos casos que suelen darse con suficiente frecuencia. 

Un flujo de pedidos se denomina flujo elemental (o flujo 
de Poisson) si el intervalo de tiempo r que transcurre entre 
dos pedidos sucesivos representa una variable aleatoria 
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distribuida en el intervalo (0, 20) con la densidad 
p (2) =0e-ox, (26) 
La loy (26) lleva también el nombre de distribución exponen- 


cial (véase la fig. 15, donde hemos representado las densi- 
dades (26) para a =1 y a = 2). 


Fig. 15 


Es fácil calcular la esperanza matemática de T: 
Mt= f xp (2) dr= f xae=ex dz, 
0 ó 
Integrando por partes (u =x y dv = ae-"*dz), obtenemos 


A, A » Calo NE] 
Mt= [—e-0x5 +] e exdi=| — > Al =>+- 
El parámetro a se denomina densidad del flujo de pedidos, 
La tórmula correspondiente al sorteo de 1 se obtiene 
fácilmente de la ecuación (23) que tiene en nuestro caso la 
forma 
1 
$ ae-adr=y. 
0 
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Calculando la integral del primer miembro, obtenemos la 
relación 
le n=y, 
de donde resulta 
1=—Lm(-y. 


Como quiera que la variable 1 — y tiene la misma dístri- 
bución que y, podemos emplear en lugar de la última fórmula 
esta otra 


1=—Liny. (27) 


5.3. Esquema de cálculo. Consideremos, pues, el fun- 
cionamiento del sistema descrito en el punto 5.1 en el caso 
de un flujo elemental de pedidos. 

Pongamos cada línea en correspondencia con una célula 
de la memoria interior de la MCE inscribiendo en la célula 
el momento en que queda libre la línea. Sea t, el momento 
en que queda libre la ¿-ésima línea. Aceptaremos que los 
cálculos se inician en el momento 7, = 0 en el que se recibe 
el primer pedido. Podemos aceptar que en este momento 
todos los t; son iguales a 7, ya que todas las líneas están 
libres. El tiempo final del cálculo será Tim =7T7, + T- 

El primer pedido va ala línea número 1. Por consiguien- 
to, esta línea permanecerá ocupada durante el tiempo tor. 
Por eso, debemos tomar en lugar de £, el valor nuevo (t;)nv = 
=T ¡+ toc, agregar uno al contador de los pedidos cum- 
plidos y pasar a considerar el segundo pedido. 

Supongamos que hemos considerado ya k pedidos. De- 
bemos entonces sortear el momento de recepción del (k + 1)- 
ésimo pedido. Con este fin escogemos el valor siguiente de y 
y calculamos, a partir de la fórmula (27), el valor siguien- 
te 1 = Ta; después determinamos el momeuto de recepción 


Tas Th + Th> 
¿Estará libre para este momento la primera línea? Para 
responder, debemos comprobar la condición 
ST ro (28) 
Si esta condición se cumple, ello significa que en el momento 


TF .+1 la línea está ya libre y puede cumplir el pedido. Debe- 
remos sustituir £, por Tay, + too, agregar uno al contador 
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de pedidos cumplidos y pasar a considerar el pedido si- 
guiente. 

Si la condición (28) no se cumple, ello significa que la 
primera línea está ocupada en el momento 7,+,. Comproba- 
mos entonces si está o no libre la segunda línea: 


to S Tar (29) 


Si se cumple la condición (29), sustituimos t, por Tr+, + 
+ toc, agregamos uno al contador de pedidos cumplidos 
y pasamos a considerar el pedido siguiente. 

Si la condición (29) no se cumple, comprobamos la condi- 
ción 


t< Tn 
Puede ocurrir que para todo i del 4 al n se tenga 


t> Tre 


o sea, puede suceder que todas las líneas estén ocupadas en 
el momento Ty. En este caso, agregamos uno al contador 
de pedidos rechazados y pasamos a considerar el pedido 
siguiente. 

Cada vez que calcnlemos Ty+,, deberemos comprobar 
además la condición de finalización del experimento 


Tr > Trino 


Si esta condición se cumple, finalizamos el experimento. 
En los contadores de pedidos cumplidos y de pedidos recha- 
zados tendremos entonces log números Heum Y Hrech 
Este experimento se repite N veces (empleando distintos 
valores de y) y se toma la media de los resultados encontrados 


N 
1 
Mitcum 27D) HMoum yr 
pan 


N 
1 
Mpireon 2 7 Y) Hrecn Gs 
je 
donde oum «p Y Brecn «y Tepresentan los valores de cum 
y de Ureen en el j-ésimo experimento. 
En la fig. 16 aparece el esquema sinóptico del programa 
correspondiente a este cálculo (si es necesario, podemos obte- 
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| Entrada de datos | 


no 
Pin del experimento Z———| — ¿Tr<Tim? = == 


a 
3 
3 


(t0nv =7h + toc >| ¿ima? 
d no 
sí 
inv =JvJ+1 
y 


Agregamos 1 al contador 
de pedidos rechazados 


y 


Agregamos 1 al contador de 
pedidos cumplidos 


Cálculo de resultados, salida 


g 


de resultados, fin del 
prota Tres =Tn + Ta 


y 


| Env=£v341 


Fig. 16 E 
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ner en la sección «Fin del experimento» los valores de Lc um «y 
y de recen y, que resultan en cada uno de los experimentos). 

5.4. Problemas más complejos. Es fácil persuadirse de 
que este mismo método se puede aplicar en sistemas mucho 
más complejos. 

Por ejemplo, to, puede ser una variable aleatoria, y no 
una magnitud constante, y tomar diferentes valores para dis- 
tintas líneas (lo que significa que las líneas se distinguen 
por sus equipos o por la calificación del personal). El es- 
quema de cálculo, en esencia, continuará siendo el mismo, 
pero habrá que sortear cada vez los valores de to¿ correspon- 
diendo a cada línea una fórmula de sorteo propia. 

También se pueden considerar los llamados sistemas 
con espera que, en lugar de rechazar inmediatamente el 
pedido que no pueden cumplir, lo conservan durante cierto 
tiempo tp (tiempo de permanencia del pedido en el sistema); 
si en el transcurso de este tiempo se libera una de las líneas, 
dicha línea atiende el pedido. 

Además se pueden considerar sistemas on los que el 
pedido en turno es atendido por la línea que se libera antes 
que otras. Por otra parte, es factible tomar en consideración 
las intermitencias aleatorias de las líneas, el tiempo aleato- 
rio de reparación de las mismas, la variación según el tiempo 
de la densidad del flujo de pedidos y otros muchos facto- 
res. 

Por supuesto, el análisis de estos sistemas exige deter- 
minados esfuerzos. Para obtener resultados de valor práctico 
debe ser ideado un modelo correcto. Esto requiere un estudio 
suficientemente minucioso de los flujos reales de pedidos, el 
cronometraje del funcionamiento de las secciones del sistema, 
ote. 

Hablando en términos generales, es preciso conocer las 
leyes probabilísticas de funcionamiento de las secciones del 
sistema. En este caso, el método de Montecarlo permite ha- 
llar las leyes probabilísticas del funcionamiento de todo el 
sistema por complejo que sea. 

Estos métodos de análisis resultan muy útiles cuando se 
trata de proyectar una empresa: en lugar de realizar un expe- 
rimento real costoso (y a veces imposible de realizar), pode- 
mos experimentar en la MCE simulando distintas varian- 
tes de organización del trabajo o del empleo de los 
equipos, 
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$06, 
ANÁLISIS 
DE LA CALIDAD 
Y DE LA SEGURIDAD 
DE PIEZAS 


6.1. Esquema clemental de análisis de la calidad. Consi- 
deremos una pieza S que consta de cierta cantidad (posible- 
mente grande) de elementos. Por ejemplo, si S es un aparato 
eléctrico, sus elementos pueden ser las resistencias (Ri), 
las capacidades (C), las lámparas, etc. Supongamos que 


22 K25%2w 


Fig. 17 


la calidad se determina por un parámetro 1 que puede ser 
calculado si se conocen los parámetros de todos los elementos 


O A CA E) 


Por ejemplo, si Ues la tensión en una sección de un circui- 
to eléctrico, podemos considerar, basándonos en la ley de 
Ohm, las ecuaciones del circuito y determinar U de estas 
ecuaciones. 

Sin embargo, los valores reales de los parámetros de log 
elementos difieren de los valores nominales. Por ejemplo, la 
resistencia representada en la fig. 17 puede tomar un valor 
cualquiera comprendido entre 20,9 y 23,1 kQ. 

Surge la cuestión siguiente: ¿cómo influyen en el valor 
de U las desviaciones que se observan entre los valores 
reales y nominales de los parámetros de todos los ele- 
mentos? 

Los límites de variación de U podrían estimarse escogien- 
do para todo elemento los valores «peores» de los parámetros. 
Pero no siempre ni mucho menos se conoce qué colección 
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de parámetros será la «peor». Además, habiendo muchos 
elementos, esta estimación resultará exagerada por cuanto 
es poco probable que todos los parámetros a la vez resulten 
«malos». 

Por eso, es más razonable tratar de estimar la esperanza 
matemática MU y la varianza DU considerando que tanto 
los parámetros de los elementos como el propio parámetro U 
son variables aleatorias. La magnitud MU representará el 
valor medio de Jen toda la partida de piezas y la magnitud 
DU indicará qué desviaciones entre Y y MU pueden pre- 
sentarse realmente. 

Recordemos que (según hemos señalado en el punto 2.2) 


MU 4f (MRa» MR» + + 5 MCay MCs + + 05 010)» 


Siendo f una función más o menos compleja, se hace impo- 
sible el cálculo analítico de la distribución de UY. A veces 
se logra determinarla experimentalmente analizando una 
gran partida de piezas elaboradas. Pero esto no siempre es 
posible y jamás puede realizarse en la fase de diseño de una 
pieza. 

Intentemos aplicar el método de Montecarlo. Para ello es 
necesario; a) conocer las características probabilísticas de 
todos los elementos y b) conocer la función f (más exacta- 
mente, saber calcular el valor de UY a partir de cualesquiera 
valores fijos de Ri, Rig + +5 Can Cima ++ 00. .de 

La distribución probabilística de los parámetros de cual- 
quier elemento se puede obtener experimentalmente anali- 
zando una gran partida de estos elementos. Dicha distri- 
bución resulta muy a menudo normal. Por eso, muchos 
investigadores parten de ello considerando, por ejemplo, 
que la resistencia del elemento representado en la fig. 17 
es una variable aleatoria normal p con la esperanza mate- 
mática Mp =22 y con 30 =1,1 (recordemos que en un 
experimento es prácticamente imposible obtener un 
valor de p que difiera de Mp en más de 30; véase la fór- 
mula (20)). 

El esquema del análisis es muy sencillo: para todo ele- 
mento se sortea primero el valor del parámetro y se calcula 
después, a partir de la fórmula (30), el valor de Y. Este expe- 
rimento se repite N veces; obtenidos los valores U,, Us, ..- 
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. +, Un, se toma aproximadamente 


N N 
DUZ > [3 UF (2 u))”]. 


Si Y toma valores grandes, podemos sustituir el factor 
ANN — 1) de la última fórmula por 1/N y entonces dicha 


¡jaja El 


Fig. 18 


fórmula resulta un simple corolario de las fórmulas (8) y (9). 
En la Estadística Matemática se demuestra que para valores 
pequeños de N conviene emplear el factor 1/(N — 1). 

6,2. Ejemplos de análisis de la seguridad. Supongamos que 
queremos estimar el tiempo medio de duración de una pieza 
aceptando que conocemos el tiempo de duración de cada uno 
de sus elementos. 

Si consideramos que el tiempo ¿, de duración de todo 
elemento es una magnitud fija, no habrá dificultad en calcu- 
lar el tiempo £ de duración de la pieza. Por ejemplo, en la 
fig. 18 viene representada esquemáticamente una pieza que 
deja de funcionar cuando falla cualquiera de sus elementos; 
en este caso tenemos 


t= mn (las tos la ta)- (31) 


En cambio para la pieza con un elemento doblado repre- 
sentada esquemáticamente en la fig. 12, tenemos 


¿= món lla; tas máx (to 5 ta); to) (32) 


ya que si falla, por ejemplo, el elemento Ni 3, la pieza con- 
tinuará funcionando a base del elemento NM 4, 
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El tiempo de duración de cualquier elemento es, en reali- 
dad, una variable aleatoria O,» Cuando decimos que el 
plazo de duración de una bombilla es de 1 000 horas, esta- 
mos hablando del valor medio MO de la variable 6; todos 
comprendemos que una bombilla se fundirá antes mientras 
que otra (de la misma partida) durará más. 

Si conocemos las densidades de distribución de 8, 
para todos los elementos que componen la pieza, podemos 


[>| a | ] 


Fig. 19 


calcular MO por el método de Montecarlo aplicando el proce- 
dimiento del punto 6.1. En efecto, para todo elemento pode- 
mos sortear el valor de la variable O»); sea tn, este valor. 
Ahora, a partir de la fórmula correspondiente (31) o (32), 
podemos hallar el valor t de la variable aleatoria O. Repi- 
tiendo este experimento un número (N) suficiente de veces, 
podemos aceptar que 


N 
1 
MO = Ñ y ty 
pan 
donde t, es el valor de £ obtenido en el j-ésimo experimento. 
Cabe subrayar que no es nada fácil determinar la distri- 
bución del plazo de duración O,,, de los elementos; por 
ejemplo, para elementos que funcionan durante largo tiempo 
se hace difícil la organización del experimento pues habrá 
que esperar a que falle una cantidad suficiente de estos 
elementos. 
6.3. Otras posibilidades del método. Los ejemplos expues- 
tos permiten ver que la metodología del análisis de la cali- 
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dad de las piezas se basa en una ide, muy simple. Es preciso 
conocer las características probubilísticas de todos los 
elementos que constituyen la pieza y saber calcular la va- 
riable que nos interesa como una función de los parámetros 
de estos elementos. La simulación permite entonces abarcar 
las propiedades aleatorias de los parámetros. 

Aparte de la esperanza matemática y de la varianza de 
la variable que nos interesa, la simulación permite obtener 
mucha más información útil. Supongamos, por ejemplo, 


E o 1 5 2 TE Y Y E 2 
AZ SA 4222= po 
1 2 3 4 5 6 7% 
Fig. 20 


que hemos obtenido una gran cantidad N de valores U,, 
Un ..., Un de la variable aleatoria Y. Basándonos en 
estos valores, podemos encontrar aproximadamente la 
densidad de distribución de Y. Esta cuestión concierne, en 
realidad, a la Estadística ya que se trata de la elaboración 
de resultados de experimentos (realizados, es verdad, en 
una MCE). Por eso, nos limitaremos a considerar un ejemplo 
concreto. 

Supongamos que hemos obtenido un total de N = 120 


valores U,, Us, » - -, Uso de la variable aleatoria U com- 
prendidos todos en los límites 
1<0,<08,5. 


Dividamos el intervalo 1 <zx<6,5 en 11 intervalos 
iguales de longitud Az = 0,5 (o en otro número de interva- 
los no muy grande ni tampoco muy pequeño) y calculemos 
cuántos valores , corresponden a cada intervalo. Hemos 
representado estas cantidades en la fig. 20. 

Dividiendo entre N = 120 estas cantidades, obtendremos 
las frecuencias correspondientes a cada intervalo. En nuestro 
caso estas frecuencias son: 0,017, 0; 0,008; 0,12; 0,20; 0,27; 
0,14; 0,16; 0,06; 0,008 y 0,017. 

En cada intervalo construimos ahora un rectángulo de 
área igual a la frecuencia de Uy correspondiente a dicho 
intervalo (fig. 21). En otras palabras, la altura de todo rectán- 
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gulo es igual a la frecuencia dividida entre Az. La figura 
escalonada así obtenida se denomina histograma. 

El histograma ofrece una aproximación de la densidad 
desconocida de la variable aleatoria Y. Por esto, el área del 
histograma comprendida entre, digamos, z = 25 yx =5,5 
nos da un valor aproximado de la probabilidad 


P (2,5< U< 5,5) = 0,95. 


El cálculo (experimento) realizado nos permite afirmar, 
por consiguiente, que con una probabilidad aproximada- 
mente igual a 0,95 la variable UY está comprendida en el 
intervalo 2,5 < U< 5,5. 

En la fig. 21 hemos representado, a título comparativo, 
la densidad de la variable aleatoria normal Z' con los pará- 
metros a == 3,85 y 6 = 0,88. Si partiendo de esta densidad 
calculamos la probabilidad de que L' pertenezca al intervalo 
2,5 < Y” < 5,5, obtendremos el valor 0,91*) bastante próxi- 
mo al valor anterior. 


1) Expliquemos cómo se efectúan los cálculos, 
Ségun (14), tenemos 
5,5 _(x-aR 


1 f zos 
e dz, 
sY 21 > 


P(2,5<Y <5,5)= 


Efectuemos en la integral la sustitución zx — a = of. Tendremos 
entonces 


donde 


Para calcular el valor de la última integral hay que recurrir a la tabla 
de la así llamada integral de probabilidades O (x) en la que vienen 
los valores de Ja función 


z Rh 
ae z 
ca po dt. 


Así encontramos 
P(2,5<U <5,5) =0,5 [0 (1,88) +0 (1,54) 0,91. 
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6.4. Observación. Lamentablemente son raros, por ahora, 
los casos en los que se aplican semejantes cálculos. Es difí- 
cil explicar esto. Posiblemente se deba a que los diseñadores 
ignoran esta posibilidad. 

Además, para poder aplicar este método al diseño de 
piezas es necesario estudiar previamente las características 


YA 


O1p 


a] 


o 1 2 3 4 5 $ 
Fig. 21 


probabilísticas de los elementos que componen la pieza. 

ste es un trabajo voluminoso. Es verdad, que, conociendo 
estas características, se puede estimar la calidad de cualquier 
pieza formada por dichos elementos. También se puede 
analizar cómo varía la calidad al sustituir un elemento por 
otro. 

Es de esperar que en los años próximos estos cálculos se 
hagan habituales. En cuanto a las características probabilís- 
ticas necesarias de los elementos, serán ofrecidas siempre por 
las empresas que los elaboren. 
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E 
ANÁLISIS DEL PASO 
DE NEUTRONES A TRAVÉS 
DE UNA PLACA 


Las leyes probabilísticas de la interacción entre una par- 
tícula elemental (neutrón, fotón, mesón, etc.) y la materia 
son conocidas. Como regla, es preciso conocer las caracterís- 
ticas macroscópicas —densidades, flujos, ote. — de procesos 
en los que participa una gran cantidad de estas partículas, 
Esta situación se asemeja a la que hemos analizado en el $ 5 
y en el $ 6 y resulta muy apropiada para la aplicación del 
método de Montecarlo. 

Posiblemente sea en la Física de neutrones donde se aplica 
con mayor frecuencia el método de Montecarlo. Analizaremos 
la variante más sencilla del problema sobre el paso de neu- 
trones a través de una placa. 

7.1. Planteamiento del problema. Supongamos que una 
placa homogénea infinita 0 < zx <h es bombardeada por 
un flujo de neutrones de energía E, siendo el ángulo de inci- 
dencia igual a 90. Al chocar con los átomos de la placa, los 
neutrones pueden esparcirse elásticamente o pueden ser 
absorbidos. Supongamos, para simplificar, que el neutrón 
conserva, al esparcirse, sn energía y que cualquier cambio 
en la dirección del neutrón, provocado por el choque con 
un átomo, es igualmente probable (esto último sucede en 
materias constituidas por átomos pesados). En la fig. 22 
hemos representado las distintas suertes que pueden correr 
los neutrones: el neutrón (a) atraviesa la placa, el neutrón 
(b) es absorbido y el neutrón (c) es reflejado. 

Es preciso calcular la probabilidad p* de que el neutrón 
atraviese la placa, la probabilidad p” de que el neutrón 
sea reflejado y la probabilidad p? de que el neutrón resulte 
absorbido. 

a interacción de los neutrones con la materia se caracte- 
riza en este caso mediante dos constantes E, y E, que se 
denominan sección de absorción y sección de dispersión, res- 
pectivamente. Los índices c y s provienen de las letras ini- 
ciales de las palabras inglesas «capture» (captura) y «scat- 
tering» (dispersión). 
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La suma de estas secciones leva el nombre de sección 
completa 


E=2.+3 


Las secciones tienen el siguiente significado físico: al 
chocar el neutrón con un átomo de la materia, la probabilidad 


Y 


Fig, 22 


de absorción es igual a E//2 y la probabilidad de esparci- 
miento es igual a Y,/2. 

El recorrido libre.) del neutrón (o sea, su recorrido entre 
un choque y otro) es una variable aleatoria que puede tomar 
cualesquiera valores positivos con la densidad de proba- 
bilidad 

p(1)= 2e-*, 


Es fácil ver que la densidad de la variable A coincide 
con la densidad (26) de la variable aleatoria t en el caso del 
flujo elemental de pedidos. Esta analogía con el punto 5.2 
nos permite escribir inmediatamente la fórmula para el 
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recorrido libre medio 


MA=1/2 
y la fórmula para el sorteo de A 
2=—41n y. 


Resta explicar cómo se escoge la dirección aleatoria del 
neutrón cuando se produce el esparcimiento. Puesto que 
existe en el problema una simetría respecto al eje x, esta 
dirección se determina plenamente con sólo indicar el ángu- 
lo p entre la dirección de la velocidad del neutrón y el eje 
Oz. Se puede demostrar*) que la condición de que todas las 
direcciones son igualmente probables es en este caso equi- 
valente a la condición de que el coseno de este ángulo p = 
= cos y esté uniformemente distribuido en el intervalo 
(—1, 1). Tomando en la fórmula (24) a =—1 y b=1, 
obtenemos la fórmula para el sorteo de pu: 

p=2y-4. 

7.2. Esquema de análisis basado en la simulación de 
trayectorias reales. Supongamos que el neutrón ha expe- 
rimentado el %-ésimo esparcimiento en el punto de abscisa 7; 
correspondiente al interior de la placa y ha comenzado a 
moverse en la dirección y. 

Sorteemos el recorrido líbre 


da = —(1/2) ln y 
y calculemos la abscisa del punto de colisión siguiente 
(fig. 23) 
Zpra = 2 + Ap 
Comprobamos la condición de que la placa ha sido atrave- 


sada 
Tr >h. 


Si esta condición se cumple, concluimos el análisis de la 
trayectoria del neutrón y agregamos uno al contador de las 
partículas que atraviesan la placa. En el caso contrario, 
comprobamos la condición de reflexión 


Lg <0. 


1) La demostración viene en el punto 9.4. 
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Si esta condición se cumple, concluimos el análisis de la 
trayectoria y agregamos uno al contador de las partículas 
reflejadas. Si tampoco se cumple esta condición, o sea, si 
0< Zr+, <h, ello significará que el neutrón experimenta 


Fig. 28 


la (k + 1)-ésima colisión en el interior de la placa; habrá 
que sortear de nuevo el «destino» del neutrón en esta colisión. 

De acuerdo con el punto 4.1 escogemos el valor siguiente 
de y y comprobamos la condición de absorción 


y< Ez. 
Si esta última desigualdad se cumple, concluimos el aná» 
lisis de la trayectoria y agregamos uno al contador de par- 
tículas absorbidas. En el caso contrario, aceptamos que 
el neutrón experimenta un esparcimiento en el punto de 
abscisa ps. Realizamos entonces el sorteo de la nueva 
dirección del neutrón 

Par =29 —1 
y volvemos a repetir todo el ciclo (tomando, claro está, 
otros valores para y). 
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Hemos escrito todos los valores de y sin índices sobren- 
tendiendo que cada valor de y se emplea una vez solamente. 
Para analizar un eslabón de la trayectoria es preciso emplear 
tres valores de y. 

Los valores iniciales para cada trayectoria son 


Zo =0 y po=1. 


Después de haber analizado de esta forma N trayectorias, 
Tesultará que N* neutrones han atravesado la placa, que N- 
neutrones han sido reflejados y que N* neutrones han quedado 
absorbidos. 

Es evidente que las probabilidades buscadas son 
aproximadamente iguales a las razones 


N+* > N- No 
Py", 1 RIP 


En la fig. 24 hemos representado el esquema sinóptico 
del programa correspondiente a este análisis, El índice ¿ 
representa el número de la trayectoria y el índice k corres- 
ponde al número de la colisión (a Jo largo de la trayectoria 
que se está analizando). 

Este método de análisis, a pesar de ser muy natural, 
tiene sus deficiencias. En particular, resulta difícil calcular 
por este método la probabilidad p* cuando ella es muy peque- 
ña. Pero precisamente con este caso nos encontramos al ana- 
lizar la protección contra la radiación. 

Existen modificaciones ingeniosas del método de Monte- 
carlo que pueden ser aplicadas también en estos casos. 

Detengámonos brevemente en la variante más sencilla 
del análisis que se basa en el empleo de los así llama- 
dos «pesos». 

7.3. Esquema de análisis basado en el empleo de pesos 
que sustituyen la absorción. Consideremos el mismo proble- 
ma sobre el paso de los neutrones. Supongamos que a lo 
largo de la trayectoria se desplaza un «tren» formado por un 
gran número w, de neutrones idénticos. A] producirse la 
Colisión en el punto de abscisa x,, un promedio de Ww, (2/2) 
de neutrones del «tren» será absorbido y un promedio de 
1w, (2,/2) de neutrones será esparcido. 

Agreguemos la cantidad w (2/2) al contador de partícu- 
las absorbidas y continuemos el análisis del «tren» esparcido 
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Entrada de datos 


sí 
Ni = Ny +1 


No =N8, +1 


y 


a T 
¿m0 ] trmtata] 
5] a | 


Cálculo de p*, p” y p? 
4 


Salida de resultados, fin 
del problema 


Fig, 26 
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aceptando que todo el «tren» restante se mueve en una misma 
dirección. 

Todas las fórmulas que hemos dado en el punto 7.2 se 
conservan. La única diferencia consiste en que, después 
de cada colisión, la cantidad de neutrones en el «tren» dis- 
minuirá, siendo 


Wrs1 04 (E4/2), 


ya que una parte del «tren», con un total de w, (2¿/2) de 
neutrones, será absorbida; además, la trayectoria jamás 
podrá terminarse con la absorción. 

La variable w, suele denominarse peso de neutrón y, en 
lugar de un «tren» formado por w, neutrones, suele hablarse 
de un neutrón de peso w,. El peso inicial w, se toma, como 
regla, igual a 1. Esto no contradice nuestra hipótesis de 
que el «tren» es grande, pues es fácil ver que wo es un factor 
común de todos los pesos w, que se obtienen al analizar una 
trayectoria. 

El esquema sinóptico del programa correspondiente a 
este análisis viene representado en la fig. 25. Como puede 
verse no resulta más complejo que el esquema de la fig. 24. 
Sin embargo, se puede demostrar!) que siempre es más con- 
veniente calcular p* por el método de este punto y no por 
el método del punto 7 

7.4, Observación. Existe una gran variedad de proce- 
dimientos de cálculo que se basan en el empleo de distintos 
pesos. No podemos detenernos aquí en estas cuestiones. 

Señalemos solamente que el método de Montecarlo es 
aplicable a problemas mucho más complejos relacionados 
con partículas elementales. Por ejemplo, al caso en que el 
medio se compone de materias distintas y tiene una estructura 
geométrica arbitraria o al caso en que la energía de la partí- 
Cula varía después de cada colisión. Así mismo, hace posible 
el análisis de otros muchos aspectos de los procesos nucleares: 
permite prever la posibilidad de que el átomo se desintegre al 
chocar con el neutrón y de que surjan neutrones nuevos; 
permito determinar el momento de inicio y la condición 
de mantenimiento de la reacción en cadena; etc. 


1) La demostración viene en el punto 9.5. 
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Entrado de datos 


A 


j=t | 
AAA y 


k=0, 20, o=1, 11) =1 


ig 
Niy=N iy +00 +] [e 
| Ape =07 (34/2) | 


Y 
Ly =NY +4 (2¿/2) | 


Jov=jvy 1 | | 


y 


7 
¿lay < N? fr =2y—1 ] 


no 4 
| rr | 


Cálculo de p*, p” y p? 


A 
Salida de resultados, ftn 
del problema 


Fig. 25 
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$8. 
CÁLCULO 
DE LA INTEGRAL DEFINIDA 


Los problemas considerados en los parágrafos 5, 6 y 7 
eran de contenido probabilístico por su esencia y la aplica- 
ción del método de Montecarlo en estos problemas parecía 
bastante natural. Ahora consideraremos un problema pura- 
mente matemático, el cálculo aproximado de la integral 
definida. 

Puesto que el cálculo de integrales definidas equivale al 
cálculo de áreas, podríamos emplear el método del punto 
1,2. Aquí explicaremos otro método, más efectivo, que 
permite construir diferentes modelos probabilísticos para 
resolver este problema por el método de Montecarlo y expli- 
caremos cómo escoger entre estos modelos los «más conve- 
nientes». 

8.1. Método de cálculo. Consideremos una función g (2) 
definida en el intervalo a <x< b. Es preciso hallar el 
valor aproximado de la integral 


b 
Tes f g(2) dz. (83) 


Tomemos una densidad de distribución py (x) arbitraria 
definida en el intervalo (a, b) (o sea, una función Py lx) 
cualquiera que cumpla las condiciones (15) y (16)). 

Además de la variable aleatoria E definida con la den- 
sidad pg (x) en el intervalo (a, hb), necesitaremos emplear 
la variable aleatoria 


n= g (Ep: (8). 
Según (18), tenemos 


b 
Mn=f le(d/p(2)1pe(2)d2=1. 


Consideremos ahora N variables aleatorias 9,, Ya, - . -, wa 
independientes e idénticas, y apliquemos a la suma de estas 
variables el teorema central del límite del punto 2.4. La 
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fórmula (21) se puede escribir en nuestro caso así 
N == 
Y D: 
P (43 m=1]<s $] =0.901. (34) 
jet 


Escojamos N valores E,, Ex, » . -, Ey. La última relación 
nos dice que siendo N suficientemente grande tendremos 


17 26 
xv 2 1 (35) 
jmi 


PE)” 


Al mismo tiempo nos dice que con una elevada probabilidad 
el error de la aproximación (35) no pasa de 


3Y Dr/N. 


8.2. Cómo escoger el esquema de cálculo. Como hemos 
visto Dar calcular la integral (33) se puede emplear una 
variable aleatoria E cualquiera definida en el intervalo (a, b). 
En todos los casos 


Mn = M lg (Ep: (8)1 = £. 
Pero la varianza Dn y, por consiguiente, también la esti- 


mación (35) del error depende de qué variable E tomemos, 
En efecto, 


» 
DI=M (1) —2=| 1e*(2)/p: (2)1d:—P. 


Se puede demostrar!) que esta expresión será mínima 
cuando pg (x=) sea proporcional a | g (2) |. 

Por supuesto, no podemos escoger fórmulas muy com- 
plejas para pz (x) ya que entonces resultará laborioso sortear 
los valores de E. Pero es posible guiarse por esta recomenda- 
ción a la hora de escoger la función p; (x) (véase el ejem- 
plo del punto 8.3). 

En la práctica el método de Montecarlo no se aplica para 
el cálculo de integrales de tipo (33) ya que existen métodos 
más exactos basados en las fórmulas de cuadratura. Pero la 
situación cambia si pasamos a considerar integrales múlti- 
ples: las fórmulas de cuadratura se hacen muy voluminosas 
mientras que ol método de Montecarlo no necesila casi 
modificaciones. 


1) La demostración viene en el punto 9.6, 
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8.3. Ejemplo numérico. Calculemos el valor aproximado de la 
integral 
m/2 
I= | senzdz. 
0 
El valor exacto de esta integral es 
m2 
f sen =dz=[—cos 2J/2=4. 
0 
Emplearemos para el cálculo dos distintas variables aleatorias E: 
una con la densidad constante p; (x) = 2/x (o sea, uniformemente 
distribuida en el intervalo (0, x/2)) y otra con la densidad lineal 
PE (2) = 8x/n2. En la fig. 26 hemos representado ambas densidades así 


8x 
y a 
1 y=senx 
2 
7 

—. 

o ca 


Flg. 26 
¡_EX-—A A a KA 
como la función integrando sen zx. La figura permite ver que la den- 
sidad lineal está en mayor consonancia con la recomendación del pun- 
to 8,2 de que es deseable que Pg (=) y sen z sean proporcionales, Por 


ello, es de esperar que el segundo procedimiento dé un resultado mejor, 
(a) Soa pg (2) = 2/x: en ol intervalo (0, 1/2). La fórmula del sorteo 


de E se puedo obtener de la fórmula (24) tomando a= 0 y b= n/2: 
= 1/2. 
La fórmula (35) da ll 


N 
xr z 
1==37 ») ni 
jul 
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Sea N = 10. Tomemos para los valores de y las ternas de cifras de 
la tabla A (multiplicándolas por 0,001). Los cálculos auxiliares apare- 
cen en la tabla 1. 

Tabla 1 


El resultado de los cálculos es 
Y == 0,952. 
(b) Sea ahora pg (2) = 8x/n”. Para sortear E empleamos la ecua- 
ción (23) 


| (82/12) dx= y; 


0 
de aquí encontramos fácilmente que 
n 
tE=z vr. 
La fórmula (35) da 
N 
eS Lt S 2ut 
EN DE" 
jal 
Sea N = 10. Los números y se escogen igual que en el caso (a). 
Los cálculos auxiliares aparecen en la tabla 2. 
Tabla 2 
¿ | 1 2 E 5 9 10 


» b.soo aso. oro. 04 os, o. 0,32 
Es |t,461 O 02.4), ese 2] 0,52 1,271 | 4,415 | 0,905 
sen ¿f 


[0, 680/0,936|0,968/0,78310,937/0,748| 0,863 | 0,751 | 0,698 | 0,868 


E 
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El resultado de los cálculos es 
1 1,016. 


El segundo procedimiento, como lo habíamos previsto, da un re- 
sultado más exacto. 


8.4. Sobre la estimación del error. En el punto 8.2 ya 
hemos señalado que el error absoluto que se obtiene al calcu- 
Jar la integral Y mediante la fórmula (35) no puede práctica- 
mente sobrepasar la magnitud 3Y Dn/N. Pero el error abso- 
luto real es, como regla, considerablemente monor. Por eso, 
para caracterizar el error suele emplearse en la práctica otra 
magnitud denominada error probable 


Spr = 0,675Y Dn/N. 


El error absoluto real depende de los números aleatorios 
que se han empleado en el cálculo y puede suceder que sea 2 
6 3 veces mayor o varias veces menor que Spr- Es decir, 8p 
indica tan sólo el orden del error y no sus límites). 


Volvamos al ejemplo del punto 8.3. Empleando los valores que 
aparecen en Jas tablas 1 y 2 y la fórmula para la varianza indicada en 
el punto 6.1, podemos calcular aproximadamente las varianzas Dn 
en ambos procedimientos ?). 

En la tabla 3 indicamos para ambos métodos los valore 
aproximados de las varianzas Dn los errores probables $pp calculados 


Tabla 3 
Método Dn Sor | dear 
(a) | 0,255 | 0.103 | 0,048 


(b) | 0,016 | 0,027 0,016 


1) Sobre el error probable véase el puuto 9,7. 
2) En el caso del método (a) tenemos 


0 10 
o 2 [3 (son EY (3 sonty)"]> 
4 et 


ea 2 
= 30 (4,604 -— 3,670) =0, 56. 
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a base de estos valores y los orrores absolutos reales 5¿¿1 que se han 


obtenido en los cálculos, 
Como vemos, 8,1 es efectivamente del mismo orden que Sor» 


APÉNDICE 


$9 
DEMOSTRACIÓN 
DE ALGUNAS 
PROPOSICIONES 


En este parágrafo insertamos la demostración de algunas 
proposiciones enunciadas en los parágrafos anteriores que o 
bien nos parecen un tanto voluminosas para el nivel que 
hemos mantenido en la exposición o bien requieren cono- 
olcnimntos más profundos de la Teoría de las probabili- 

lades. 

9.1, Argumentación del método de Neumann para el 
sorteo de una variable aleatoria (punto 4.3). El punto alea- 
torio T' está uniformemente distribuido en el rectángulo abcd 
(fig. 27) cuyo área es igual a M, (b — a)!). La probabilidad 
de que el punto P aparezca por debajo de la curva y = 
=p (2) y no sea rechazado de la consideración, es igual a 
la razón entre las áreas 


b 
NES 
Mb) Mba 


La probabilidad de que el punto T' aparezca por debajo 
de la curva y = p (x) y corresponda al intervalo a << 


En el caso del método (b) tenemos 


10 10 
mn senEyy2 4 sen Ey y2 
y (3 y)- 
4 
5% 

1) Véase el punto 9.3. 


(6,875 —6,777) 0,016. 
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< b' es también igual a la razón entre las áreas 


Por consiguiente, entre todos los valores que hemos esco- 
gido para E, la parte de los valores que pertenecen al inter- 


YA 
A EAS: TEA A 
y=p(x) 
o o ar b” e 


Pig. 27 
— 
valo (a', b”) es igual al cociente 
pr 
$ p(x) dz 1 
a A ez 
anta =P 
e 
que es lo que se quería demostrar. 
9.2. Densidad de distribución de la variable ' =a + 
+ 0% (punto 4.4). Se supone que la variable £ es una varia- 
ble normal con la esperanza matemática Mi =60 y la 
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varianza Dí = 1 de modo que su densidad es 


ss ESE, 
Prix) Vas 
Para calcular la densidad de distribución de la variable 
T', escojamos dos números arbitrarios 1, < x, y calcule- 
mos la probabilidad 


P(x<U <zj=P [,<a4- 0d < ze) 


=P [ as 
Por consiguiente, 
(a-ayo _% 
Pl<E<a)= 337 e * dz. 
o/o 


Realizando en la última integral la sustitución a” =a + 0%, 
obtenemos 


P<U<mp= 


1 
V2xo 
de donde resulta (véase (14)) que la variable ¿' es una varia- 
ble normal con los parámetros Mi' =a y DE' =0?. 

9.3. Distribución uniforme de los puntos (y”, y”) en 
el cuadrado (punto 4.5). Puesto que las coordenadas del 
punto (y”, y””) son independientes, su densidad p (zx, y) 
es igual al producto de las densidades 


P (x, y) = Py (2) Pyr (y). 
Cada una de estas densidades es igual idénticamente a 1. 
Por consiguiente, p (zx, y) =1 (para 0O<zx<1 y 0< 
< y < 1). Pero esto significa precisamente que el punto 
(y, y””) está uniformemente distribuido en el cuadrado de 
dimensión 1%. 

9.4. Selección de la dirección aleatoria (punto 7.1). La 
dirección se puede determinar mediante un vector de lon- 
gitud 1 que arranca del origen de coordenadas. Los extre- 
mos de estos vectores se encuentran en la superficie esférica 
de radio 1. Las palabras «cualquier dirección es igualmente 
probable» significan que el extremo del vector representa 
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un punto aleatorio 2 uniformemente distribuido sobre la 
superficie esférica; es decir, la probabilidad de que Q per- 


A 2. a ds 
tenezca a cualquier elemento d$S de la superficie es igual a Za 


Fig. 28 


Tomemos en la superficie esférica las coordenadas esféóri- 
cas (p, 1p) con el eje polar Oz (fig. 28). Tendremos entonces 


dS = sen y dy dy (36) 
con la particularidad de que (O<py<xm, 0O< p< 2n. 


Sea p (p, y) la densidad del punto aleatorio (q, Y). 
Debemos tener 


Pp (9, Y) de dy = dS/án; 


de aquí y de la relación (36) se deduce que 


== (87) 


plo) =+ 


7 


A partir de la densidad conjunta de q y Y, podemos calcular 
fácilmente las densidades de estas variables: 


La 

p9=| po pap, (38) 
o 

n(y)= j PO, Y do=37. (39) 


o 


La igualdad p (p, Y) = Po (9) Po (4) demuestra que y y 
y son independientes. 

Es obvio que y está uniformemente distribuida en el 
intervalo (0, 21) y, por eso, la fórmula para el sorteo de p es 


Y = 21y. (40) 


La fórmula para el sorteo de q se obtiene de la ecua- 
ción (23): 


e 
+ senado=, 
0 


de donde resulta 
cos q =1 — 2y. (41) 


Las fórmulas (40) y (41) permiten escoger (sortear) la di- 
rección aleatoria. Por supuesto, los valores de y en estas 
fórmulas deben ser distintos. 

La última fórmula del punto 7.1 se obtiene de la fórmula 
(41) si se toma 1 — y en lugar de y; pero estas variables tienen 
la misma distribución. 

9.5. Ventaja del método de análisis basado en el empleo 
de pesos (punlo 7.3). Consideremos las variables aleatorias 
v y v' iguales a la cantidad (peso) de neutrones de una misma 
trayectoria que atraviesan la placa, con la particularidad 
de que la primera variable corresponde al método del punto 
7.2 y la segunda, al método del punto 7.3. 

Según el contenido físico del problema tenemos 


Mv = Mv” = p*. 


(Véase Ja demostración rigurosa de esta afirmación en el 
libro [3].) 
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Puesto que v puede tomar sólo dos valores, O y 1, la dis- 
tribución de v es 
4 0 
a O 
Tomando en cuenta que v? = v, es fácil ver que Dv = 


= p* — (p*)*. Es evidente que la variable y” puede tomar, 
aparte del valor 0, un conjunto infinito de valores wy = 1, 


Wi, Waz, . + .y Uh, » . . Por eso, su distribución es 
Li “ep 
lo € Ya -»- Y ---Q 


Los valores concretos de q,, como veremos, no nos harán falta, 
Cualesquiera que sean podemos representar la varianza en 
la forma 


DY = 2 wkga —(p*). 
h=0 


Observando que todos los w, < 1 y que 2Q10ndn =M= 


== p*, obtenemos la desigualdad Dv” < p* — (p*)? = Dv. 

El hecho de que la varianza de v' sea siempre menor que 
la varianza de v significa que, al calcular p*, el método del 
punto 7.3 dará siempre mejores resultados que el método 
del punto 7.2. La misma conclusión podemos hacer en cuanto 
al cálculo de p- y también al cáclulo de p” si la absorción 
no es muy grande. 


9.6. Selección óptima de E (punto 8.2). En el punto 8.2 
hemos obtenido la expresión para la varianza Dn. Para 
hallar el mínimo de esta expresión, siendo arbitraria la selec- 
ción de p¿(z), recurriremos a la siguiente desigualdad 
conocida del Análisis 


Ñ luto 122] <f ve (2) de- Í 2% (2) de. 


Tomando u =8 (I/V p; (2) y v=VP:(. obtenemos de 
esta desigualdad 


ñ let122]<) La dz | mtajar= LElbaz 
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Por consiguiente, 
2 
D>[ 100142]. (42) 


Resta demostrar que la cota inferior se alcanza cuando 
Ps (x) es proporcional a | g (1) |. 
Sea 


pa(x) = LEC... (43) 
$18 (2)1 dz 


Es fácil calcular que para esta densidad py (x) se tiene 
Y a 
a “2 
Fi coimt1da=[ | 16 (01 az] 
2 A 


de modo que la varianza Dn resulta efectivamente igual al 
segundo miembro de (42). 

Obsorvemos que es prácticamente imposible emplear en el 
cálculo la densidad «más convenientes (43): deberíamos 


conocer para ello el valor de la integral j lg (2) | dz; pero el 


a 

cálculo de esta integral es, de hecho, un problema equiva- 

lente al problema que tenemos planteado acerca del cálculo de 
b 


la integral f g (x) dz. Por eso, nos hemos limitado a dar la 


a 
recomendación indicada en el punto 8.2. 

9.7. Definición del error probable (punto 8.4). Sea £ la 
variable aleatoria normal definida en el punto 2.3. Es fácil 
calcular que para r = 0,6750 y cualesquiera a y o se tiene 


ajr 
j Pi (2) dx= 0,5. 
a=r 
De aquí se deduce que 
P lb—a|<r) =P ([8-2a]>r) =055, 
o sea, que las desviaciones mayores que r y las desviaciones 
menores que 7 son igualmente probables. Por oso, la mag- 
nitud r se denomina error probable de la variable aleatoria É. 
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En el punto 8.1 hemos realizado-el cálculo aproximado de 
la variable normal p = (1/N) (n, +12 +... + Ny). Re- 
sultó que su esperanza matemática es a = Mp = Í y que 
su varianza es 0? = Dp = Dn/N. Por eso, el error probable 
de la variable p es aproximadamente igual a 0,675/ Dn/N. 


$ 10. 
SOBRE LOS NÚMEROS 
SEUDOALEATORIOS 


Los algoritmos que se emplean para obtener números 
seudoaleatorios tienen, en su mayoría, la forma 


Var =P (11)- (44) 


Si se ha dado el número inicial po, todos los números siguien- 
tes Yi, Va » +. Se calculan aplicando la fórmula (44) 
con k=0, 1, 2, ... Ambos algoritmos considerados en el 
punto 3.3 son también de la forma (44) sólo, en lugar de 
señalar la forma analítica de la función y = F (x), se ha 
indicado el conjunto de operaciones que debe efectuarse con 
el argumento zx para obtener el valor de y. 

10.1. Cómo debe escogerse la función 7 (x). El ejemplo 
que sigue explica en qué consiste una de las dificultades 
principales que se plantea al escoger F (x). 


Ejemplo. Demostremos que la función y = F' (z) cuyo gráfico 
aa en la fig. 29 no se puede emplear para obtener, a partir de la 
fórmula (44), números seudoaleatoríos, 

En efecto, consideremos en el cuadrado 0. <r<1,0<y<1 
los puntos de coordenadas 

We 72d» (Vas Ya, (Vs, VO)» + > ++ 
donde y, = F (1), Pz F (va), ve = F (yo). + - - Todos estos puntos 
ge encuentran sobre la curva y = F (x) y esto, por supuesto, está muy 
mal ya que los puntos aleatorios verdaderos deben cubrir uniforme- 
mente todo el cuadrado (véase el punto 9.3). 


Este ejemplo pormite ver que el empleo de la función 
y =PF(z) on la fórmula (44) puede resultar exitoso sólo 
si el gráfico de esta función recubre con bastante densidad 
todo el cuadrado. 

Este es el caso, por ejemplo, de la función 


y = (ez), (45) 
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donde g es un número muy grande y (z) representa la parte 
fraccionaria del número z. En la fig. 30 aparece el gráfico 


Fig. 29 


de la función (45) correspondiente a g = 20. El lector podrá 
hacerse una idea de lo que representará este gráfico si g = 5%, 

10.2. Método de congruencias (método de residuos). El 
algoritmo que más se emplea para obtener números seudo- 
aleatorios fue propuesto por D. Lehmer. Se basa en la fun- 
ción (45), pero suele construirse de una manera algo distinta 
que permite realizarlo fácilmente en las MCE. 

Se determina una sucesión de números enteros mz toman- 
do my, =1 y calculando todos los números sucesivos my, 


May + 2 partir de una misma fórmula 
Mass =5"m, (mod 200), (46) 
donde k=0, 1, 2, ... Basándose en los números m;, se 


calculan los números seudoaleatorios 
Y=2 my. (47) 


7 


La fórmula (46) significa que el número my+, es igual al 
resto que se obtiene al dividir 51m, entre 2%, lón la teoría 
de las congruencias (véase cualquier texto de Teoría de los 


y=[20x) 


Fig. 30 


números) este resto lleva el nombre de residuo positivo mó- 
dulo 2 mínimo. Por eso, el algoritmo suele denominarse 
método de congruencias o método de residuos. 

Es muy fácil realizar las fórmulas (46) y (47) en las 
MCE que funcionan con números de 40 órdenes aplicando 
la operación de multiplicación con el número duplicado 
de órdenes!): hay que emplear las cifras correspondientes a 
los órdenes menores del producto. 

El período de la sucesión ma es igual a 23, 

En el caso de las MCE que funcionan con número de 36 órdones, se 
toman en las fórmulas (46) y (47) los números 515 y 2% en lugar de 51? 
y 210; el período de la sucesión my es igual entonces a 2%, 


1) La MCE Strelá (véase el punto 3.3) no puede 
realizar esta operación. 


52452 42499 
76773 97528 
04825 82134 
87113 84718 
57616 

0,2005 1,4922 
1,1609 —0,6690 
0,5864 —0,9245 
0,4425 0,2863 
0,9516 4,7708 
0,5863 0,8574 
1,1572 0,9990 
—0,4428 — —0,5564 
0,3924 7984 
0,8319 0,4270 
0,9780 —0,7679 


Tabla A. 400 cifras aleatorias D 


66 455 66 434 56 558 12 332 
42 502 99 224 88 955 53 758 
85 181 38 967 33 181 72.664 
88059 89 342 67248 09 082 
89 688 78416 27589 99 528 
33 346 83 935 79 130 90410 
27256 66 447 25 731 37525 
80 317 75120 45 904 75604 
45863 24 520 19976 04 925 
38 132 64 294 15218 49 286 
Tabla B. 88 valores normales 2 

—0,0077 0,0348 1,0423 —1,8149 
1,5893 0,5816 1,8818 0,7390 
,0904 1,5068  —1,4147 0,2776 
4,2809 4043 0,6379  —0,4428 
2,8854 0,4686 1,4664 4,6852 
—0,5557 0,8115  —0,2676  —1,2496 
—0,1032 0,5405  —0,6022 0,0093 
-—0,5098  —1,1929 —0,0572  —0,5061 
0,6141 —41,3598 1,4943  —0,4406 
—0,8888 0,4167 — -—0,8513 1,1054 
0,8960 0,5454 — —0,7165 0,8563 


94 317 
91 644 
53 807 
12311 
14 480 
45 420 
16 287 
70 492 
07 824 
89571 


1,4803 
—0,2736 
0,4012 
—2,3006 
—0,9690 
—41,2125 
0,2119 
0,4557 
0,2033 
1,2287 
4,1630 


Tablas 


1) Las cifras aleatorías imitan los valores de la variable aleatoria con la distribución (22) (véase el punto 3,1). 


3) Los valores normales imitan los valores de la variable alealoría normai (gaussiana) Y con los parámetros a = 0 


yo=í. 
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